src/MEDCoupling/MEDCoupling1GTUMesh.cxx

   1 // Copyright (C) 2007-2013  CEA/DEN, EDF R&D
   2 //
   3 // This library is free software; you can redistribute it and/or
   4 // modify it under the terms of the GNU Lesser General Public
   5 // License as published by the Free Software Foundation; either
   6 // version 2.1 of the License.
   7 //
   8 // This library is distributed in the hope that it will be useful,
   9 // but WITHOUT ANY WARRANTY; without even the implied warranty of
  10 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  11 // Lesser General Public License for more details.
  12 //
  13 // You should have received a copy of the GNU Lesser General Public
  14 // License along with this library; if not, write to the Free Software
  15 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
  16 //
  17 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
  18 //
  19 // Author : Anthony Geay (CEA/DEN)
  20
  21 #include "MEDCoupling1GTUMesh.hxx"
  22 #include "MEDCouplingUMesh.hxx"
  23 #include "MEDCouplingFieldDouble.hxx"
  24
  25 #include "SplitterTetra.hxx"
  26
  27 using namespace ParaMEDMEM;
  28
  29 MEDCoupling1GTUMesh::MEDCoupling1GTUMesh()
  30 {
  31 }
  32
  33 MEDCoupling1GTUMesh::MEDCoupling1GTUMesh(const char *name, const INTERP_KERNEL::CellModel& cm):_cm(&cm)
  34 {
  35   setName(name);
  36 }
  37
  38 MEDCoupling1GTUMesh::MEDCoupling1GTUMesh(const MEDCoupling1GTUMesh& other, bool recDeepCpy):MEDCouplingPointSet(other,recDeepCpy),_cm(other._cm)
  39 {
  40 }
  41
  42 MEDCoupling1GTUMesh *MEDCoupling1GTUMesh::New(const char *name, INTERP_KERNEL::NormalizedCellType type) throw(INTERP_KERNEL::Exception)
  43 {
  44   if(type==INTERP_KERNEL::NORM_ERROR)
  45     throw INTERP_KERNEL::Exception("MEDCoupling1GTUMesh::New : NORM_ERROR is not a valid type to be used as base geometric type for a mesh !");
  46   const INTERP_KERNEL::CellModel& cm=INTERP_KERNEL::CellModel::GetCellModel(type);
  47   if(!cm.isDynamic())
  48     return MEDCoupling1SGTUMesh::New(name,type);
  49   else
  50     return MEDCoupling1DGTUMesh::New(name,type);
  51 }
  52
  53 MEDCoupling1GTUMesh *MEDCoupling1GTUMesh::New(const MEDCouplingUMesh *m) throw(INTERP_KERNEL::Exception)
  54 {
  55   if(!m)
  56     throw INTERP_KERNEL::Exception("MEDCoupling1GTUMesh::New : input mesh is null !");
  57   std::set<INTERP_KERNEL::NormalizedCellType> gts(m->getAllGeoTypes());
  58   if(gts.size()!=1)
  59     throw INTERP_KERNEL::Exception("MEDCoupling1GTUMesh::New : input mesh must have exactly one geometric type !");
  60   const INTERP_KERNEL::CellModel& cm=INTERP_KERNEL::CellModel::GetCellModel(*gts.begin());
  61   if(!cm.isDynamic())
  62     return MEDCoupling1SGTUMesh::New(m);
  63   else
  64     return MEDCoupling1DGTUMesh::New(m);
  65 }
  66
  67 const INTERP_KERNEL::CellModel& MEDCoupling1GTUMesh::getCellModel() const throw(INTERP_KERNEL::Exception)
  68 {
  69   return *_cm;
  70 }
  71
  72 INTERP_KERNEL::NormalizedCellType MEDCoupling1GTUMesh::getCellModelEnum() const throw(INTERP_KERNEL::Exception)
  73 {
  74   return _cm->getEnum();
  75 }
  76
  77 int MEDCoupling1GTUMesh::getMeshDimension() const
  78 {
  79   return (int)_cm->getDimension();
  80 }
  81
  82 /*!
  83  * This method returns a newly allocated array containing cell ids (ascendingly sorted) whose geometric type are equal to type.
  84  * This method does not throw exception if geometric type \a type is not in \a this.
  85  * This method throws an INTERP_KERNEL::Exception if meshdimension of \b this is not equal to those of \b type.
  86  * The coordinates array is not considered here.
  87  *
  88  * \param [in] type the geometric type
  89  * \return cell ids in this having geometric type \a type.
  90  */
  91 DataArrayInt *MEDCoupling1GTUMesh::giveCellsWithType(INTERP_KERNEL::NormalizedCellType type) const throw(INTERP_KERNEL::Exception)
  92 {
  93   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret=DataArrayInt::New();
  94   if(type==getCellModelEnum())
  95     ret->alloc(getNumberOfCells(),1);
  96   else
  97     ret->alloc(0,1);
  98   ret->iota();
  99   return ret.retn();
 100 }
 101
 102 /*!
 103  * Returns nb of cells having the geometric type \a type. No throw if no cells in \a this has the geometric type \a type.
 104  */
 105 int MEDCoupling1GTUMesh::getNumberOfCellsWithType(INTERP_KERNEL::NormalizedCellType type) const
 106 {
 107   return type==getCellModelEnum()?getNumberOfCells():0;
 108 }
 109
 110 /*!
 111  * Returns a type of a cell by its id.
 112  *  \param [in] cellId - the id of the cell of interest.
 113  *  \return INTERP_KERNEL::NormalizedCellType - enumeration item describing the cell type.
 114  *  \throw If \a cellId is invalid. Valid range is [0, \a this->getNumberOfCells() ).
 115  */
 116 INTERP_KERNEL::NormalizedCellType MEDCoupling1GTUMesh::getTypeOfCell(int cellId) const
 117 {
 118   if(cellId>=0 && cellId<getNumberOfCells())
 119     return getCellModelEnum();
 120   std::ostringstream oss; oss << "MEDCoupling1GTUMesh::getTypeOfCell : Requesting type of cell #" << cellId << " but it should be in [0," << getNumberOfCells() << ") !";
 121   throw INTERP_KERNEL::Exception(oss.str().c_str());
 122 }
 123
 124 /*!
 125  * Returns a set of all cell types available in \a this mesh.
 126  * \return std::set<INTERP_KERNEL::NormalizedCellType> - the set of cell types.
 127  * \warning this method does not throw any exception even if \a this is not defined.
 128  */
 129 std::set<INTERP_KERNEL::NormalizedCellType> MEDCoupling1GTUMesh::getAllGeoTypes() const
 130 {
 131   std::set<INTERP_KERNEL::NormalizedCellType> ret;
 132   ret.insert(getCellModelEnum());
 133   return ret;
 134 }
 135
 136 /*!
 137  * This method expects that \a this is sorted by types. If not an exception will be thrown.
 138  * This method returns in the same format as code (see MEDCouplingUMesh::checkTypeConsistencyAndContig or MEDCouplingUMesh::splitProfilePerType) how
 139  * \a this is composed in cell types.
 140  * The returned array is of size 3*n where n is the number of different types present in \a this.
 141  * For every k in [0,n] ret[3*k+2]==-1 because it has no sense here.
 142  * This parameter is kept only for compatibility with other methode listed above.
 143  */
 144 std::vector<int> MEDCoupling1GTUMesh::getDistributionOfTypes() const throw(INTERP_KERNEL::Exception)
 145 {
 146   std::vector<int> ret(3);
 147   ret[0]=(int)getCellModelEnum(); ret[1]=getNumberOfCells(); ret[2]=-1;
 148   return ret;
 149 }
 150
 151 /*!
 152  * This method is the opposite of MEDCouplingUMesh::checkTypeConsistencyAndContig method. Given a list of cells in \a profile it returns a list of sub-profiles sorted by geo type.
 153  * The result is put in the array \a idsPerType. In the returned parameter \a code, foreach i \a code[3*i+2] refers (if different from -1) to a location into the \a idsPerType.
 154  * This method has 1 input \a profile and 3 outputs \a code \a idsInPflPerType and \a idsPerType.
 155  *
 156  * \param [out] code is a vector of size 3*n where n is the number of different geometric type in \a this \b reduced to the profile \a profile. \a code has exactly the same semantic than in MEDCouplingUMesh::checkTypeConsistencyAndContig method.
 157  * \param [out] idsInPflPerType is a vector of size of different geometric type in the subpart defined by \a profile of \a this ( equal to \a code.size()/3). For each i,
 158  *              \a idsInPflPerType[i] stores the tuple ids in \a profile that correspond to the geometric type code[3*i+0]
 159  * \param [out] idsPerType is a vector of size of different sub profiles needed to be defined to represent the profile \a profile for a given geometric type.
 160  *              This vector can be empty in case of all geometric type cells are fully covered in ascending in the given input \a profile.
 161  *
 162  * \warning for performance reasons no deep copy will be performed, if \a profile can been used as this in output parameters \a idsInPflPerType and \a idsPerType.
 163  *
 164  * \throw if \a profile has not exactly one component. It throws too, if \a profile contains some values not in [0,getNumberOfCells()) or if \a this is not fully defined
 165  *
 166  *  \b Example1: <br>
 167  *          - Before \a this has 3 cells \a profile contains [0,1,2]
 168  *          - After \a code contains [NORM_...,nbCells,-1], \a idsInPflPerType [[0,1,2]] and \a idsPerType is empty <br>
 169  *
 170  *  \b Example2: <br>
 171  *          - Before \a this has 3 cells \a profile contains [1,2]
 172  *          - After \a code contains [NORM_...,nbCells,0], \a idsInPflPerType [[0,1]] and \a idsPerType is [[1,2]] <br>
 173
 174  */
 175 void MEDCoupling1GTUMesh::splitProfilePerType(const DataArrayInt *profile, std::vector<int>& code, std::vector<DataArrayInt *>& idsInPflPerType, std::vector<DataArrayInt *>& idsPerType) const throw(INTERP_KERNEL::Exception)
 176 {
 177   if(!profile)
 178     throw INTERP_KERNEL::Exception("MEDCoupling1GTUMesh::splitProfilePerType : input profile is NULL !");
 179   if(profile->getNumberOfComponents()!=1)
 180     throw INTERP_KERNEL::Exception("MEDCoupling1GTUMesh::splitProfilePerType : input profile should have exactly one component !");
 181   int nbTuples=profile->getNumberOfTuples();
 182   int nbOfCells=getNumberOfCells();
 183   code.resize(3); idsInPflPerType.resize(1);
 184   code[0]=(int)getCellModelEnum(); code[1]=nbTuples;
 185   idsInPflPerType.resize(1);
 186   if(profile->isIdentity() && nbTuples==nbOfCells)
 187     {
 188       code[2]=-1;
 189       idsInPflPerType[0]=const_cast<DataArrayInt *>(profile); idsInPflPerType[0]->incrRef();
 190       idsPerType.clear();
 191       return ;
 192     }
 193   code[2]=0;
 194   profile->checkAllIdsInRange(0,nbOfCells);
 195   idsPerType.resize(1);
 196   idsPerType[0]=const_cast<DataArrayInt *>(profile); idsPerType[0]->incrRef();
 197   idsInPflPerType[0]=DataArrayInt::Range(0,nbTuples,1);
 198 }
 199
 200 /*!
 201  * This method tries to minimize at most the number of deep copy.
 202  * So if \a idsPerType is not empty it can be returned directly (without copy, but with ref count incremented) in return.
 203  *
 204  * \sa MEDCouplingUMesh::checkTypeConsistencyAndContig
 205  */
 206 DataArrayInt *MEDCoupling1GTUMesh::checkTypeConsistencyAndContig(const std::vector<int>& code, const std::vector<const DataArrayInt *>& idsPerType) const throw(INTERP_KERNEL::Exception)
 207 {
 208   int nbOfCells=getNumberOfCells();
 209   if(code.size()!=3)
 210     throw INTERP_KERNEL::Exception("MEDCoupling1GTUMesh::checkTypeConsistencyAndContig : invalid input code should be exactly of size 3 !");
 211   if(code[0]!=(int)getCellModelEnum())
 212     {
 213       std::ostringstream oss; oss << "MEDCoupling1GTUMesh::checkTypeConsistencyAndContig : Mismatch of geometric type ! Asking for " << code[0] << " whereas the geometric type is \a this is " << getCellModelEnum() << " (" << _cm->getRepr() << ") !";
 214       throw INTERP_KERNEL::Exception(oss.str().c_str());
 215     }
 216   if(code[2]==-1)
 217     {
 218       if(code[1]==nbOfCells)
 219         return 0;
 220       else
 221         {
 222           std::ostringstream oss; oss << "MEDCoupling1GTUMesh::checkTypeConsistencyAndContig : mismatch between the number of cells in this (" << nbOfCells << ") and the number of non profile (" << code[1] << ") !";
 223           throw INTERP_KERNEL::Exception(oss.str().c_str());
 224         }
 225     }
 226   if(code[2]!=0)
 227     throw INTERP_KERNEL::Exception("MEDCoupling1GTUMesh::checkTypeConsistencyAndContig : single geo type mesh ! 0 or -1 is expected at pos #2 of input code !");
 228   if(idsPerType.size()!=1)
 229     throw INTERP_KERNEL::Exception("MEDCoupling1GTUMesh::checkTypeConsistencyAndContig : input code points to DataArrayInt #0 whereas the size of idsPerType is not equal to 1 !");
 230   const DataArrayInt *pfl=idsPerType[0];
 231   if(!pfl)
 232     throw INTERP_KERNEL::Exception("MEDCoupling1GTUMesh::checkTypeConsistencyAndContig : the input code points to a NULL DataArrayInt at rank 0 !");
 233   if(pfl->getNumberOfComponents()!=1)
 234     throw INTERP_KERNEL::Exception("MEDCoupling1GTUMesh::checkTypeConsistencyAndContig : input profile should have exactly one component !");
 235   pfl->checkAllIdsInRange(0,nbOfCells);
 236   pfl->incrRef();
 237   return const_cast<DataArrayInt *>(pfl);
 238 }
 239
 240 void MEDCoupling1GTUMesh::writeVTKLL(std::ostream& ofs, const std::string& cellData, const std::string& pointData, DataArrayByte *byteData) const throw(INTERP_KERNEL::Exception)
 241 {
 242   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> m=buildUnstructured();
 243   m->writeVTKLL(ofs,cellData,pointData,byteData);
 244 }
 245
 246 std::string MEDCoupling1GTUMesh::getVTKDataSetType() const throw(INTERP_KERNEL::Exception)
 247 {
 248   return std::string("UnstructuredGrid");
 249 }
 250
 251 std::size_t MEDCoupling1GTUMesh::getHeapMemorySize() const
 252 {
 253   return MEDCouplingPointSet::getHeapMemorySize();
 254 }
 255
 256 bool MEDCoupling1GTUMesh::isEqualIfNotWhy(const MEDCouplingMesh *other, double prec, std::string& reason) const throw(INTERP_KERNEL::Exception)
 257 {
 258   if(!MEDCouplingPointSet::isEqualIfNotWhy(other,prec,reason))
 259     return false;
 260   if(!other)
 261     throw INTERP_KERNEL::Exception("MEDCoupling1GTUMesh::isEqualIfNotWhy : input other pointer is null !");
 262   const MEDCoupling1GTUMesh *otherC=dynamic_cast<const MEDCoupling1GTUMesh *>(other);
 263   if(!otherC)
 264     {
 265       reason="mesh given in input is not castable in MEDCouplingSGTUMesh !";
 266       return false;
 267     }
 268   if(_cm!=otherC->_cm)
 269     {
 270       reason="mismatch in geometric type !";
 271       return false;
 272     }
 273   return true;
 274 }
 275
 276 bool MEDCoupling1GTUMesh::isEqualWithoutConsideringStr(const MEDCouplingMesh *other, double prec) const
 277 {
 278   if(!MEDCouplingPointSet::isEqualWithoutConsideringStr(other,prec))
 279     return false;
 280   if(!other)
 281     throw INTERP_KERNEL::Exception("MEDCoupling1GTUMesh::isEqualWithoutConsideringStr : input other pointer is null !");
 282   const MEDCoupling1GTUMesh *otherC=dynamic_cast<const MEDCoupling1GTUMesh *>(other);
 283   if(!otherC)
 284     return false;
 285   if(_cm!=otherC->_cm)
 286     return false;
 287   return true;
 288 }
 289
 290 void MEDCoupling1GTUMesh::checkCoherency() const throw(INTERP_KERNEL::Exception)
 291 {
 292   MEDCouplingPointSet::checkCoherency();
 293 }
 294
 295 DataArrayDouble *MEDCoupling1GTUMesh::getBarycenterAndOwner() const
 296 {
 297   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> m=buildUnstructured();
 298   MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> ret=m->getBarycenterAndOwner();
 299   return ret.retn();
 300 }
 301
 302 MEDCouplingFieldDouble *MEDCoupling1GTUMesh::getMeasureField(bool isAbs) const
 303 {
 304   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> m=buildUnstructured();
 305   MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> ret=m->getMeasureField(isAbs);
 306   ret->setMesh(this);
 307   return ret.retn();
 308 }
 309
 310 MEDCouplingFieldDouble *MEDCoupling1GTUMesh::getMeasureFieldOnNode(bool isAbs) const
 311 {
 312   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> m=buildUnstructured();
 313   MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> ret=m->getMeasureFieldOnNode(isAbs);
 314   ret->setMesh(this);
 315   return ret.retn();
 316 }
 317
 318 /*!
 319  * to improve perf !
 320  */
 321 int MEDCoupling1GTUMesh::getCellContainingPoint(const double *pos, double eps) const
 322 {
 323   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> m=buildUnstructured();
 324   return m->getCellContainingPoint(pos,eps);
 325 }
 326
 327 MEDCouplingFieldDouble *MEDCoupling1GTUMesh::buildOrthogonalField() const
 328 {
 329   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> m=buildUnstructured();
 330   MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> ret=m->buildOrthogonalField();
 331   ret->setMesh(this);
 332   return ret.retn();
 333 }
 334
 335 DataArrayInt *MEDCoupling1GTUMesh::getCellsInBoundingBox(const double *bbox, double eps) const
 336 {
 337   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> m=buildUnstructured();
 338   return m->getCellsInBoundingBox(bbox,eps);
 339 }
 340
 341 DataArrayInt *MEDCoupling1GTUMesh::getCellsInBoundingBox(const INTERP_KERNEL::DirectedBoundingBox& bbox, double eps)
 342 {
 343   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> m=buildUnstructured();
 344   return m->getCellsInBoundingBox(bbox,eps);
 345 }
 346
 347 MEDCouplingPointSet *MEDCoupling1GTUMesh::buildFacePartOfMySelfNode(const int *start, const int *end, bool fullyIn) const
 348 {
 349   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> m=buildUnstructured();
 350   return m->buildFacePartOfMySelfNode(start,end,fullyIn);
 351 }
 352
 353 DataArrayInt *MEDCoupling1GTUMesh::findBoundaryNodes() const
 354 {
 355   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> m=buildUnstructured();
 356   return m->findBoundaryNodes();
 357 }
 358
 359 MEDCouplingPointSet *MEDCoupling1GTUMesh::buildBoundaryMesh(bool keepCoords) const
 360 {
 361   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> m=buildUnstructured();
 362   return m->buildBoundaryMesh(keepCoords);
 363 }
 364
 365 void MEDCoupling1GTUMesh::findCommonCells(int compType, int startCellId, DataArrayInt *& commonCellsArr, DataArrayInt *& commonCellsIArr) const throw(INTERP_KERNEL::Exception)
 366 {
 367   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> m=buildUnstructured();
 368   m->findCommonCells(compType,startCellId,commonCellsArr,commonCellsIArr);
 369 }
 370
 371 int MEDCoupling1GTUMesh::getNodalConnectivityLength() const throw(INTERP_KERNEL::Exception)
 372 {
 373   const DataArrayInt *c1(getNodalConnectivity());
 374   if(!c1)
 375     throw INTERP_KERNEL::Exception("MEDCoupling1GTUMesh::getNodalConnectivityLength : no connectivity set !");
 376   if(c1->getNumberOfComponents()!=1)
 377     throw INTERP_KERNEL::Exception("MEDCoupling1GTUMesh::getNodalConnectivityLength : Nodal connectivity array set must have exactly one component !");
 378   if(!c1->isAllocated())
 379     throw INTERP_KERNEL::Exception("MEDCoupling1GTUMesh::getNodalConnectivityLength : Nodal connectivity array must be allocated !");
 380   return c1->getNumberOfTuples();
 381 }
 382
 383 /*!
 384  * This method aggregates all the meshes in \a parts to put them in a single unstructured mesh (those returned).
 385  * The order of cells is the returned instance is those in the order of instances in \a parts.
 386  *
 387  * \param [in] parts - all not null parts of single geo type meshes to be aggreagated having the same mesh dimension and same coordinates.
 388  * \return MEDCouplingUMesh * - new object to be dealt by the caller.
 389  *
 390  * \throw If one element is null in \a parts.
 391  * \throw If not all the parts do not have the same mesh dimension.
 392  * \throw If not all the parts do not share the same coordinates.
 393  * \throw If not all the parts have their connectivity set properly.
 394  * \throw If \a parts is empty.
 395  */
 396 MEDCouplingUMesh *MEDCoupling1GTUMesh::AggregateOnSameCoordsToUMesh(const std::vector< const MEDCoupling1GTUMesh *>& parts) throw(INTERP_KERNEL::Exception)
 397 {
 398   if(parts.empty())
 399     throw INTERP_KERNEL::Exception("MEDCoupling1GTUMesh::AggregateOnSameCoordsToUMesh : input parts vector is empty !");
 400   const MEDCoupling1GTUMesh *firstPart(parts[0]);
 401   if(!firstPart)
 402     throw INTERP_KERNEL::Exception("MEDCoupling1GTUMesh::AggregateOnSameCoordsToUMesh : the first instance in input parts is null !");
 403   const DataArrayDouble *coords(firstPart->getCoords());
 404   int meshDim(firstPart->getMeshDimension());
 405   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> ret(MEDCouplingUMesh::New(firstPart->getName().c_str(),meshDim)); ret->setDescription(firstPart->getDescription().c_str());
 406   ret->setCoords(coords);
 407   int nbOfCells(0),connSize(0);
 408   for(std::vector< const MEDCoupling1GTUMesh *>::const_iterator it=parts.begin();it!=parts.end();it++)
 409     {
 410       if(!(*it))
 411         throw INTERP_KERNEL::Exception("MEDCoupling1GTUMesh::AggregateOnSameCoordsToUMesh : presence of null pointer in input vector !");
 412       if((*it)->getMeshDimension()!=meshDim)
 413         throw INTERP_KERNEL::Exception("MEDCoupling1GTUMesh::AggregateOnSameCoordsToUMesh : all the instances in input vector must have same mesh dimension !");
 414       if((*it)->getCoords()!=coords)
 415         throw INTERP_KERNEL::Exception("MEDCoupling1GTUMesh::AggregateOnSameCoordsToUMesh : all the instances must share the same coordinates pointer !");
 416       nbOfCells+=(*it)->getNumberOfCells();
 417       connSize+=(*it)->getNodalConnectivityLength();
 418     }
 419   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> conn(DataArrayInt::New()),connI(DataArrayInt::New());
 420   connI->alloc(nbOfCells+1,1); conn->alloc(connSize+nbOfCells,1);
 421   int *c(conn->getPointer()),*ci(connI->getPointer()); *ci=0;
 422   for(std::vector< const MEDCoupling1GTUMesh *>::const_iterator it=parts.begin();it!=parts.end();it++)
 423     {
 424       int curNbCells((*it)->getNumberOfCells());
 425       int geoType((int)(*it)->getCellModelEnum());
 426       const int *cinPtr((*it)->getNodalConnectivity()->begin());
 427       const MEDCoupling1SGTUMesh *ps(dynamic_cast<const MEDCoupling1SGTUMesh *>(*it));
 428       const MEDCoupling1DGTUMesh *pd(dynamic_cast<const MEDCoupling1DGTUMesh *>(*it));
 429       if(ps && !pd)
 430         {
 431           int nNodesPerCell(ps->getNumberOfNodesPerCell());
 432           for(int i=0;i<curNbCells;i++,ci++,cinPtr+=nNodesPerCell)
 433             {
 434               *c++=geoType;
 435               c=std::copy(cinPtr,cinPtr+nNodesPerCell,c);
 436               ci[1]=ci[0]+nNodesPerCell+1;
 437             }
 438         }
 439       else if(!ps && pd)
 440         {
 441           const int *ciinPtr(pd->getNodalConnectivityIndex()->begin());
 442           for(int i=0;i<curNbCells;i++,ci++,ciinPtr++)
 443             {
 444               *c++=geoType;
 445               c=std::copy(cinPtr+ciinPtr[0],cinPtr+ciinPtr[1],c);
 446               ci[1]=ci[0]+ciinPtr[1]-ciinPtr[0]+1;
 447             }
 448         }
 449       else
 450         throw INTERP_KERNEL::Exception("MEDCoupling1GTUMesh::AggregateOnSameCoordsToUMesh : presence of instance which type is not in [MEDCoupling1SGTUMesh,MEDCoupling1DGTUMesh] !");
 451     }
 452   ret->setConnectivity(conn,connI,true);
 453   return ret.retn();
 454 }
 455
 456 //==
 457
 458 MEDCoupling1SGTUMesh::MEDCoupling1SGTUMesh(const MEDCoupling1SGTUMesh& other, bool recDeepCpy):MEDCoupling1GTUMesh(other,recDeepCpy),_conn(other._conn)
 459 {
 460   if(recDeepCpy)
 461     {
 462       const DataArrayInt *c(other._conn);
 463       if(c)
 464         _conn=c->deepCpy();
 465     }
 466 }
 467
 468 MEDCoupling1SGTUMesh::MEDCoupling1SGTUMesh(const char *name, const INTERP_KERNEL::CellModel& cm):MEDCoupling1GTUMesh(name,cm)
 469 {
 470 }
 471
 472 MEDCoupling1SGTUMesh::MEDCoupling1SGTUMesh()
 473 {
 474 }
 475
 476 MEDCoupling1SGTUMesh *MEDCoupling1SGTUMesh::New()
 477 {
 478   return new MEDCoupling1SGTUMesh;
 479 }
 480
 481 MEDCoupling1SGTUMesh *MEDCoupling1SGTUMesh::New(const char *name, INTERP_KERNEL::NormalizedCellType type) throw(INTERP_KERNEL::Exception)
 482 {
 483   if(type==INTERP_KERNEL::NORM_ERROR)
 484     throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::New : NORM_ERROR is not a valid type to be used as base geometric type for a mesh !");
 485   const INTERP_KERNEL::CellModel& cm=INTERP_KERNEL::CellModel::GetCellModel(type);
 486   if(cm.isDynamic())
 487     {
 488       std::ostringstream oss; oss << "MEDCoupling1SGTUMesh::New : the input geometric type " << cm.getRepr() << " is dynamic ! Only static types are allowed here !";
 489       throw INTERP_KERNEL::Exception(oss.str().c_str());
 490     }
 491   return new MEDCoupling1SGTUMesh(name,cm);
 492 }
 493
 494 MEDCoupling1SGTUMesh *MEDCoupling1SGTUMesh::New(const MEDCouplingUMesh *m) throw(INTERP_KERNEL::Exception)
 495 {
 496   if(!m)
 497     throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::New : input mesh is null !");
 498   std::set<INTERP_KERNEL::NormalizedCellType> gts(m->getAllGeoTypes());
 499   if(gts.size()!=1)
 500     throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::New : input mesh must have exactly one geometric type !");
 501   int geoType((int)*gts.begin());
 502   MEDCouplingAutoRefCountObjectPtr<MEDCoupling1SGTUMesh> ret(MEDCoupling1SGTUMesh::New(m->getName().c_str(),*gts.begin()));
 503   ret->setCoords(m->getCoords()); ret->setDescription(m->getDescription().c_str());
 504   int nbCells(m->getNumberOfCells());
 505   int nbOfNodesPerCell(ret->getNumberOfNodesPerCell());
 506   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> conn(DataArrayInt::New()); conn->alloc(nbCells*nbOfNodesPerCell,1);
 507   int *c(conn->getPointer());
 508   const int *cin(m->getNodalConnectivity()->begin()),*ciin(m->getNodalConnectivityIndex()->begin());
 509   for(int i=0;i<nbCells;i++,ciin++)
 510     {
 511       if(cin[ciin[0]]==geoType)
 512         {
 513           if(ciin[1]-ciin[0]==nbOfNodesPerCell+1)
 514             c=std::copy(cin+ciin[0]+1,cin+ciin[1],c);
 515           else
 516             {
 517               std::ostringstream oss; oss << "MEDCoupling1SGTUMesh::New(const MEDCouplingUMesh *m) : something is wrong in the input mesh at cell #" << i << " ! The size of cell is not those expected (" << nbOfNodesPerCell << ") !";
 518               throw INTERP_KERNEL::Exception(oss.str().c_str());
 519             }
 520         }
 521       else
 522         {
 523           std::ostringstream oss; oss << "MEDCoupling1SGTUMesh::New(const MEDCouplingUMesh *m) : something is wrong in the input mesh at cell #" << i << " ! The geometric type is not those expected !";
 524           throw INTERP_KERNEL::Exception(oss.str().c_str());
 525         }
 526     }
 527   ret->setNodalConnectivity(conn);
 528   return ret.retn();
 529 }
 530
 531 MEDCoupling1SGTUMesh *MEDCoupling1SGTUMesh::clone(bool recDeepCpy) const
 532 {
 533   return new MEDCoupling1SGTUMesh(*this,recDeepCpy);
 534 }
 535
 536 /*!
 537  * This method behaves mostly like MEDCoupling1SGTUMesh::deepCpy method, except that only nodal connectivity arrays are deeply copied.
 538  * The coordinates are shared between \a this and the returned instance.
 539  *
 540  * \return MEDCouplingUMesh * - A new object instance holding the copy of \a this (deep for connectivity, shallow for coordiantes)
 541  * \sa MEDCoupling1SGTUMesh::deepCpy
 542  */
 543 MEDCouplingPointSet *MEDCoupling1SGTUMesh::deepCpyConnectivityOnly() const throw(INTERP_KERNEL::Exception)
 544 {
 545   checkCoherency();
 546   MEDCouplingAutoRefCountObjectPtr<MEDCoupling1SGTUMesh> ret(clone(false));
 547   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> c(_conn->deepCpy());
 548   ret->setNodalConnectivity(c);
 549   return ret.retn();
 550 }
 551
 552 void MEDCoupling1SGTUMesh::shallowCopyConnectivityFrom(const MEDCouplingPointSet *other) throw(INTERP_KERNEL::Exception)
 553 {
 554   if(!other)
 555     throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::shallowCopyConnectivityFrom : input pointer is null !");
 556   const MEDCoupling1SGTUMesh *otherC=dynamic_cast<const MEDCoupling1SGTUMesh *>(other);
 557   if(!otherC)
 558     throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::shallowCopyConnectivityFrom : input pointer is not an MEDCoupling1SGTUMesh instance !");
 559   setNodalConnectivity(otherC->getNodalConnectivity());
 560 }
 561
 562 void MEDCoupling1SGTUMesh::updateTime() const
 563 {
 564   MEDCoupling1GTUMesh::updateTime();
 565   const DataArrayInt *c(_conn);
 566   if(c)
 567     updateTimeWith(*c);
 568 }
 569
 570 std::size_t MEDCoupling1SGTUMesh::getHeapMemorySize() const
 571 {
 572   std::size_t ret=0;
 573   const DataArrayInt *c(_conn);
 574   if(c)
 575     ret+=c->getHeapMemorySize();
 576   return MEDCoupling1GTUMesh::getHeapMemorySize()+ret;
 577 }
 578
 579 MEDCouplingMesh *MEDCoupling1SGTUMesh::deepCpy() const
 580 {
 581   return clone(true);
 582 }
 583
 584 bool MEDCoupling1SGTUMesh::isEqualIfNotWhy(const MEDCouplingMesh *other, double prec, std::string& reason) const throw(INTERP_KERNEL::Exception)
 585 {
 586   if(!other)
 587     throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::isEqualIfNotWhy : input other pointer is null !");
 588   std::ostringstream oss; oss.precision(15);
 589   const MEDCoupling1SGTUMesh *otherC=dynamic_cast<const MEDCoupling1SGTUMesh *>(other);
 590   if(!otherC)
 591     {
 592       reason="mesh given in input is not castable in MEDCoupling1SGTUMesh !";
 593       return false;
 594     }
 595   if(!MEDCoupling1GTUMesh::isEqualIfNotWhy(other,prec,reason))
 596     return false;
 597   const DataArrayInt *c1(_conn),*c2(otherC->_conn);
 598   if(c1==c2)
 599     return true;
 600   if(!c1 || !c2)
 601     {
 602       reason="in connectivity of single static geometric type exactly one among this and other is null !";
 603       return false;
 604     }
 605   if(!c1->isEqualIfNotWhy(*c2,reason))
 606     {
 607       reason.insert(0,"Nodal connectivity DataArrayInt differ : ");
 608       return false;
 609     }
 610   return true;
 611 }
 612
 613 bool MEDCoupling1SGTUMesh::isEqualWithoutConsideringStr(const MEDCouplingMesh *other, double prec) const
 614 {
 615   if(!other)
 616     throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::isEqualWithoutConsideringStr : input other pointer is null !");
 617   const MEDCoupling1SGTUMesh *otherC=dynamic_cast<const MEDCoupling1SGTUMesh *>(other);
 618   if(!otherC)
 619     return false;
 620   if(!MEDCoupling1GTUMesh::isEqualWithoutConsideringStr(other,prec))
 621     return false;
 622   const DataArrayInt *c1(_conn),*c2(otherC->_conn);
 623   if(c1==c2)
 624     return true;
 625   if(!c1 || !c2)
 626     return false;
 627   if(!c1->isEqualWithoutConsideringStr(*c2))
 628     return false;
 629   return true;
 630 }
 631
 632 void MEDCoupling1SGTUMesh::checkCoherencyOfConnectivity() const throw(INTERP_KERNEL::Exception)
 633 {
 634   const DataArrayInt *c1(_conn);
 635   if(c1)
 636     {
 637       if(c1->getNumberOfComponents()!=1)
 638         throw INTERP_KERNEL::Exception("Nodal connectivity array is expected to be with number of components set to one !");
 639       if(c1->getInfoOnComponent(0)!="")
 640         throw INTERP_KERNEL::Exception("Nodal connectivity array is expected to have no info on its single component !");
 641       c1->checkAllocated();
 642     }
 643   else
 644     throw INTERP_KERNEL::Exception("Nodal connectivity array not defined !");
 645 }
 646
 647 void MEDCoupling1SGTUMesh::checkCoherency() const throw(INTERP_KERNEL::Exception)
 648 {
 649   MEDCouplingPointSet::checkCoherency();
 650   checkCoherencyOfConnectivity();
 651 }
 652
 653 void MEDCoupling1SGTUMesh::checkCoherency1(double eps) const throw(INTERP_KERNEL::Exception)
 654 {
 655   checkCoherency();
 656   const DataArrayInt *c1(_conn);
 657   int nbOfTuples=c1->getNumberOfTuples();
 658   int nbOfNodesPerCell=(int)_cm->getNumberOfNodes();
 659   if(nbOfTuples%nbOfNodesPerCell!=0)
 660     {
 661       std::ostringstream oss; oss << "MEDCoupling1SGTUMesh::checkCoherency1 : the nb of tuples in conn is " << nbOfTuples << " and number of nodes per cell is " << nbOfNodesPerCell << ". But " << nbOfTuples << "%" << nbOfNodesPerCell << " !=0 !";
 662       throw INTERP_KERNEL::Exception(oss.str().c_str());
 663     }
 664   int nbOfNodes=getNumberOfNodes();
 665   int nbOfCells=nbOfTuples/nbOfNodesPerCell;
 666   const int *w(c1->begin());
 667   for(int i=0;i<nbOfCells;i++)
 668     for(int j=0;j<nbOfNodesPerCell;j++,w++)
 669       {
 670         if(*w<0 || *w>=nbOfNodes)
 671           {
 672             std::ostringstream oss; oss << "At node #" << j << " of  cell #" << i << ", is equal to " << *w << " must be in [0," << nbOfNodes << ") !";
 673             throw INTERP_KERNEL::Exception(oss.str().c_str());
 674           }
 675       }
 676 }
 677
 678 void MEDCoupling1SGTUMesh::checkCoherency2(double eps) const throw(INTERP_KERNEL::Exception)
 679 {
 680   checkCoherency1(eps);
 681 }
 682
 683 int MEDCoupling1SGTUMesh::getNumberOfCells() const
 684 {
 685   int nbOfTuples=getNodalConnectivityLength();
 686   int nbOfNodesPerCell=getNumberOfNodesPerCell();
 687   if(nbOfTuples%nbOfNodesPerCell!=0)
 688     {
 689       std::ostringstream oss; oss << "MEDCoupling1SGTUMesh:getNumberOfCells: : the nb of tuples in conn is " << nbOfTuples << " and number of nodes per cell is " << nbOfNodesPerCell << ". But " << nbOfTuples << "%" << nbOfNodesPerCell << " !=0 !";
 690       throw INTERP_KERNEL::Exception(oss.str().c_str());
 691     }
 692   return nbOfTuples/nbOfNodesPerCell;
 693 }
 694
 695 int MEDCoupling1SGTUMesh::getNumberOfNodesInCell(int cellId) const throw(INTERP_KERNEL::Exception)
 696 {
 697   return getNumberOfNodesPerCell();
 698 }
 699
 700 int MEDCoupling1SGTUMesh::getNumberOfNodesPerCell() const throw(INTERP_KERNEL::Exception)
 701 {
 702   checkNonDynamicGeoType();
 703   return (int)_cm->getNumberOfNodes();
 704 }
 705
 706 DataArrayInt *MEDCoupling1SGTUMesh::computeNbOfNodesPerCell() const throw(INTERP_KERNEL::Exception)
 707 {
 708   checkNonDynamicGeoType();
 709   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret=DataArrayInt::New();
 710   ret->alloc(getNumberOfCells(),1);
 711   ret->fillWithValue((int)_cm->getNumberOfNodes());
 712   return ret.retn();
 713 }
 714
 715 DataArrayInt *MEDCoupling1SGTUMesh::computeNbOfFacesPerCell() const throw(INTERP_KERNEL::Exception)
 716 {
 717   checkNonDynamicGeoType();
 718   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret=DataArrayInt::New();
 719   ret->alloc(getNumberOfCells(),1);
 720   ret->fillWithValue((int)_cm->getNumberOfSons());
 721   return ret.retn();
 722 }
 723
 724 DataArrayInt *MEDCoupling1SGTUMesh::computeEffectiveNbOfNodesPerCell() const throw(INTERP_KERNEL::Exception)
 725 {
 726   checkNonDynamicGeoType();
 727   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret=DataArrayInt::New();
 728   int nbCells(getNumberOfCells());
 729   ret->alloc(nbCells,1);
 730   int *retPtr(ret->getPointer());
 731   int nbNodesPerCell(getNumberOfNodesPerCell());
 732   const int *conn(_conn->begin());
 733   for(int i=0;i<nbCells;i++,conn+=nbNodesPerCell,retPtr++)
 734     {
 735       std::set<int> s(conn,conn+nbNodesPerCell);
 736       *retPtr=(int)s.size();
 737     }
 738   return ret.retn();
 739 }
 740
 741 void MEDCoupling1SGTUMesh::getNodeIdsOfCell(int cellId, std::vector<int>& conn) const
 742 {
 743   int sz=getNumberOfNodesPerCell();
 744   conn.resize(sz);
 745   if(cellId>=0 && cellId<getNumberOfCells())
 746     std::copy(_conn->begin()+cellId*sz,_conn->begin()+(cellId+1)*sz,conn.begin());
 747   else
 748     {
 749       std::ostringstream oss; oss << "MEDCoupling1SGTUMesh::getNodeIdsOfCell : request for cellId #" << cellId << " must be in [0," << getNumberOfCells() << ") !";
 750       throw INTERP_KERNEL::Exception(oss.str().c_str());
 751     }
 752 }
 753
 754 void MEDCoupling1SGTUMesh::checkNonDynamicGeoType() const throw(INTERP_KERNEL::Exception)
 755 {
 756   if(_cm->isDynamic())
 757     throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::checkNonDynamicGeoType : internal error ! the internal geo type is dynamic ! should be static !");
 758 }
 759
 760 std::string MEDCoupling1SGTUMesh::simpleRepr() const
 761 {
 762   static const char msg0[]="No coordinates specified !";
 763   std::ostringstream ret;
 764   ret << "Single static geometic type (" << _cm->getRepr() << ") unstructured mesh with name : \"" << getName() << "\"\n";
 765   ret << "Description of mesh : \"" << getDescription() << "\"\n";
 766   int tmpp1,tmpp2;
 767   double tt=getTime(tmpp1,tmpp2);
 768   ret << "Time attached to the mesh [unit] : " << tt << " [" << getTimeUnit() << "]\n";
 769   ret << "Iteration : " << tmpp1  << " Order : " << tmpp2 << "\n";
 770   ret << "Mesh dimension : " << getMeshDimension() << "\nSpace dimension : ";
 771   if(_coords!=0)
 772     {
 773       const int spaceDim=getSpaceDimension();
 774       ret << spaceDim << "\nInfo attached on space dimension : ";
 775       for(int i=0;i<spaceDim;i++)
 776         ret << "\"" << _coords->getInfoOnComponent(i) << "\" ";
 777       ret << "\n";
 778     }
 779   else
 780     ret << msg0 << "\n";
 781   ret << "Number of nodes : ";
 782   if(_coords!=0)
 783     ret << getNumberOfNodes() << "\n";
 784   else
 785     ret << msg0 << "\n";
 786   ret << "Number of cells : ";
 787   if((const DataArrayInt *)_conn)
 788     {
 789       if(_conn->isAllocated())
 790         {
 791           if(_conn->getNumberOfComponents()==1)
 792             ret << getNumberOfCells() << "\n";
 793           else
 794             ret << "Nodal connectivity array specified and allocated but with not exactly one component !" << "\n";
 795         }
 796       else
 797         ret << "Nodal connectivity array specified but not allocated !" << "\n";
 798     }
 799   else
 800     ret << "No connectivity specified !" << "\n";
 801   ret << "Cell type : " << _cm->getRepr() << "\n";
 802   return ret.str();
 803 }
 804
 805 std::string MEDCoupling1SGTUMesh::advancedRepr() const
 806 {
 807   std::ostringstream ret;
 808   ret << simpleRepr();
 809   ret << "\nCoordinates array : \n___________________\n\n";
 810   if(_coords)
 811     _coords->reprWithoutNameStream(ret);
 812   else
 813     ret << "No array set !\n";
 814   ret << "\n\nConnectivity array : \n____________________\n\n";
 815   //
 816   if((const DataArrayInt *)_conn)
 817     {
 818       if(_conn->isAllocated())
 819         {
 820           if(_conn->getNumberOfComponents()==1)
 821             {
 822              int nbOfCells=getNumberOfCells();
 823              int sz=getNumberOfNodesPerCell();
 824              const int *connPtr=_conn->begin();
 825              for(int i=0;i<nbOfCells;i++,connPtr+=sz)
 826                {
 827                  ret << "Cell #" << i << " : ";
 828                  std::copy(connPtr,connPtr+sz,std::ostream_iterator<int>(ret," "));
 829                  ret << "\n";
 830                }
 831             }
 832           else
 833             ret << "Nodal connectivity array specified and allocated but with not exactly one component !" << "\n";
 834         }
 835       else
 836         ret << "Nodal connectivity array specified but not allocated !" << "\n";
 837     }
 838   else
 839     ret << "No connectivity specified !" << "\n";
 840   return ret.str();
 841 }
 842
 843 DataArrayDouble *MEDCoupling1SGTUMesh::computeIsoBarycenterOfNodesPerCell() const throw(INTERP_KERNEL::Exception)
 844 {
 845   MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> ret=DataArrayDouble::New();
 846   int spaceDim=getSpaceDimension();
 847   int nbOfCells=getNumberOfCells();//checkCoherency()
 848   int nbOfNodes=getNumberOfNodes();
 849   ret->alloc(nbOfCells,spaceDim);
 850   double *ptToFill=ret->getPointer();
 851   const double *coor=_coords->begin();
 852   const int *nodal=_conn->begin();
 853   int sz=getNumberOfNodesPerCell();
 854   double coeff=1./(double)sz;
 855   for(int i=0;i<nbOfCells;i++,ptToFill+=spaceDim)
 856     {
 857       std::fill(ptToFill,ptToFill+spaceDim,0.);
 858       for(int j=0;j<sz;j++,nodal++)
 859         if(*nodal>=0 && *nodal<nbOfNodes)
 860           std::transform(coor+spaceDim*nodal[0],coor+spaceDim*(nodal[0]+1),ptToFill,ptToFill,std::plus<double>());
 861         else
 862           {
 863             std::ostringstream oss; oss << "MEDCoupling1SGTUMesh::computeIsoBarycenterOfNodesPerCell : on cell #" << i << " presence of nodeId #" << *nodal << " should be in [0," <<   nbOfNodes << ") !";
 864             throw INTERP_KERNEL::Exception(oss.str().c_str());
 865           }
 866       std::transform(ptToFill,ptToFill+spaceDim,ptToFill,std::bind2nd(std::multiplies<double>(),coeff));
 867     }
 868   return ret.retn();
 869 }
 870
 871 void MEDCoupling1SGTUMesh::renumberCells(const int *old2NewBg, bool check) throw(INTERP_KERNEL::Exception)
 872 {
 873   int nbCells=getNumberOfCells();
 874   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> o2n=DataArrayInt::New();
 875   o2n->useArray(old2NewBg,false,C_DEALLOC,nbCells,1);
 876   if(check)
 877     o2n=o2n->checkAndPreparePermutation();
 878   //
 879   const int *conn=_conn->begin();
 880   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> n2o=o2n->invertArrayO2N2N2O(nbCells);
 881   const int *n2oPtr=n2o->begin();
 882   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> newConn=DataArrayInt::New();
 883   newConn->alloc(_conn->getNumberOfTuples(),1);
 884   newConn->copyStringInfoFrom(*_conn);
 885   int sz=getNumberOfNodesPerCell();
 886   //
 887   int *newC=newConn->getPointer();
 888   for(int i=0;i<nbCells;i++,newC+=sz)
 889     {
 890       int pos=n2oPtr[i];
 891       std::copy(conn+pos*sz,conn+(pos+1)*sz,newC);
 892     }
 893   _conn=newConn;
 894 }
 895
 896 /*!
 897  * Keeps from \a this only cells which constituing point id are in the ids specified by [\a begin,\a end).
 898  * The resulting cell ids are stored at the end of the 'cellIdsKept' parameter.
 899  * Parameter \a fullyIn specifies if a cell that has part of its nodes in ids array is kept or not.
 900  * If \a fullyIn is true only cells whose ids are \b fully contained in [\a begin,\a end) tab will be kept.
 901  *
 902  * \param [in] begin input start of array of node ids.
 903  * \param [in] end input end of array of node ids.
 904  * \param [in] fullyIn input that specifies if all node ids must be in [\a begin,\a end) array to consider cell to be in.
 905  * \param [in,out] cellIdsKeptArr array where all candidate cell ids are put at the end.
 906  */
 907 void MEDCoupling1SGTUMesh::fillCellIdsToKeepFromNodeIds(const int *begin, const int *end, bool fullyIn, DataArrayInt *&cellIdsKeptArr) const
 908 {
 909   int nbOfCells=getNumberOfCells();
 910   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> cellIdsKept=DataArrayInt::New(); cellIdsKept->alloc(0,1);
 911   int tmp=-1;
 912   int sz=_conn->getMaxValue(tmp); sz=std::max(sz,0)+1;
 913   std::vector<bool> fastFinder(sz,false);
 914   for(const int *work=begin;work!=end;work++)
 915     if(*work>=0 && *work<sz)
 916       fastFinder[*work]=true;
 917   const int *conn=_conn->begin();
 918   int nbNodesPerCell=getNumberOfNodesPerCell();
 919   for(int i=0;i<nbOfCells;i++,conn+=nbNodesPerCell)
 920     {
 921       int ref=0,nbOfHit=0;
 922       for(int j=0;j<nbNodesPerCell;j++)
 923         if(conn[j]>=0)
 924           {
 925             ref++;
 926             if(fastFinder[conn[j]])
 927               nbOfHit++;
 928           }
 929       if((ref==nbOfHit && fullyIn) || (nbOfHit!=0 && !fullyIn))
 930         cellIdsKept->pushBackSilent(i);
 931     }
 932   cellIdsKeptArr=cellIdsKept.retn();
 933 }
 934
 935 MEDCouplingMesh *MEDCoupling1SGTUMesh::mergeMyselfWith(const MEDCouplingMesh *other) const
 936 {
 937   if(other->getType()!=SINGLE_STATIC_GEO_TYPE_UNSTRUCTURED)
 938     throw INTERP_KERNEL::Exception("Merge of umesh only available with umesh single static geo type each other !");
 939   const MEDCoupling1SGTUMesh *otherC=static_cast<const MEDCoupling1SGTUMesh *>(other);
 940   return Merge1SGTUMeshes(this,otherC);
 941 }
 942
 943 MEDCouplingUMesh *MEDCoupling1SGTUMesh::buildUnstructured() const throw(INTERP_KERNEL::Exception)
 944 {
 945   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> ret=MEDCouplingUMesh::New(getName().c_str(),getMeshDimension());
 946   ret->setCoords(getCoords());
 947   const int *nodalConn=_conn->begin();
 948   int nbCells=getNumberOfCells();
 949   int nbNodesPerCell=getNumberOfNodesPerCell();
 950   int geoType=(int)getCellModelEnum();
 951   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> c=DataArrayInt::New(); c->alloc(nbCells*(nbNodesPerCell+1),1);
 952   int *cPtr=c->getPointer();
 953   for(int i=0;i<nbCells;i++,nodalConn+=nbNodesPerCell)
 954     {
 955       *cPtr++=geoType;
 956       cPtr=std::copy(nodalConn,nodalConn+nbNodesPerCell,cPtr);
 957     }
 958   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> cI=DataArrayInt::Range(0,(nbCells+1)*(nbNodesPerCell+1),nbNodesPerCell+1);
 959   ret->setConnectivity(c,cI,true);
 960   return ret.retn();
 961 }
 962
 963 DataArrayInt *MEDCoupling1SGTUMesh::simplexize(int policy) throw(INTERP_KERNEL::Exception)
 964 {
 965   switch(policy)
 966     {
 967     case 0:
 968       return simplexizePol0();
 969     case 1:
 970       return simplexizePol1();
 971     case (int) INTERP_KERNEL::PLANAR_FACE_5:
 972       return simplexizePlanarFace5();
 973     case (int) INTERP_KERNEL::PLANAR_FACE_6:
 974       return simplexizePlanarFace6();
 975     default:
 976       throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::simplexize : unrecognized policy ! Must be :\n  - 0 or 1 (only available for meshdim=2) \n  - PLANAR_FACE_5, PLANAR_FACE_6  (only for meshdim=3)");
 977     }
 978 }
 979
 980 /// @cond INTERNAL
 981
 982 struct MEDCouplingAccVisit
 983 {
 984   MEDCouplingAccVisit():_new_nb_of_nodes(0) { }
 985   int operator()(int val) { if(val!=-1) return _new_nb_of_nodes++; else return -1; }
 986   int _new_nb_of_nodes;
 987 };
 988
 989 /// @endcond
 990
 991 /*!
 992  * Finds nodes not used in any cell and returns an array giving a new id to every node
 993  * by excluding the unused nodes, for which the array holds -1. The result array is
 994  * a mapping in "Old to New" mode.
 995  *  \param [out] nbrOfNodesInUse - number of node ids present in the nodal connectivity.
 996  *  \return DataArrayInt * - a new instance of DataArrayInt. Its length is \a
 997  *          this->getNumberOfNodes(). It holds for each node of \a this mesh either -1
 998  *          if the node is unused or a new id else. The caller is to delete this
 999  *          array using decrRef() as it is no more needed.
1000  *  \throw If the coordinates array is not set.
1001  *  \throw If the nodal connectivity of cells is not defined.
1002  *  \throw If the nodal connectivity includes an invalid id.
1003  */
1004 DataArrayInt *MEDCoupling1SGTUMesh::getNodeIdsInUse(int& nbrOfNodesInUse) const throw(INTERP_KERNEL::Exception)
1005 {
1006   nbrOfNodesInUse=-1;
1007   int nbOfNodes=getNumberOfNodes();
1008   int nbOfCells=getNumberOfCells();
1009   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret=DataArrayInt::New();
1010   ret->alloc(nbOfNodes,1);
1011   int *traducer=ret->getPointer();
1012   std::fill(traducer,traducer+nbOfNodes,-1);
1013   const int *conn=_conn->begin();
1014   int nbNodesPerCell=getNumberOfNodesPerCell();
1015   for(int i=0;i<nbOfCells;i++)
1016     for(int j=0;j<nbNodesPerCell;j++,conn++)
1017       if(*conn>=0 && *conn<nbOfNodes)
1018         traducer[*conn]=1;
1019       else
1020         {
1021           std::ostringstream oss; oss << "MEDCoupling1SGTUMesh::getNodeIdsInUse : In cell #" << i  << " presence of node id " <<  conn[j] << " not in [0," << nbOfNodes << ") !";
1022           throw INTERP_KERNEL::Exception(oss.str().c_str());
1023         }
1024   nbrOfNodesInUse=(int)std::count(traducer,traducer+nbOfNodes,1);
1025   std::transform(traducer,traducer+nbOfNodes,traducer,MEDCouplingAccVisit());
1026   return ret.retn();
1027 }
1028
1029 /*!
1030  * Changes ids of nodes within the nodal connectivity arrays according to a permutation
1031  * array in "Old to New" mode. The node coordinates array is \b not changed by this method.
1032  * This method is a generalization of shiftNodeNumbersInConn().
1033  *  \warning This method performs no check of validity of new ids. **Use it with care !**
1034  *  \param [in] newNodeNumbersO2N - a permutation array, of length \a
1035  *         this->getNumberOfNodes(), in "Old to New" mode.
1036  *         See \ref MEDCouplingArrayRenumbering for more info on renumbering modes.
1037  *  \throw If the nodal connectivity of cells is not defined.
1038  */
1039 void MEDCoupling1SGTUMesh::renumberNodesInConn(const int *newNodeNumbersO2N)
1040 {
1041   getNumberOfCells();//only to check that all is well defined.
1042   _conn->transformWithIndArr(newNodeNumbersO2N,newNodeNumbersO2N+getNumberOfNodes());
1043   updateTime();
1044 }
1045
1046 MEDCoupling1SGTUMesh *MEDCoupling1SGTUMesh::Merge1SGTUMeshes(const MEDCoupling1SGTUMesh *mesh1, const MEDCoupling1SGTUMesh *mesh2) throw(INTERP_KERNEL::Exception)
1047 {
1048   std::vector<const MEDCoupling1SGTUMesh *> tmp(2);
1049   tmp[0]=const_cast<MEDCoupling1SGTUMesh *>(mesh1); tmp[1]=const_cast<MEDCoupling1SGTUMesh *>(mesh2);
1050   return Merge1SGTUMeshes(tmp);
1051 }
1052
1053 MEDCoupling1SGTUMesh *MEDCoupling1SGTUMesh::Merge1SGTUMeshes(std::vector<const MEDCoupling1SGTUMesh *>& a) throw(INTERP_KERNEL::Exception)
1054 {
1055   std::size_t sz=a.size();
1056   if(sz==0)
1057     return Merge1SGTUMeshesLL(a);
1058   for(std::size_t ii=0;ii<sz;ii++)
1059     if(!a[ii])
1060       {
1061         std::ostringstream oss; oss << "MEDCoupling1SGTUMesh::Merge1SGTUMeshes : item #" << ii << " in input array of size "<< sz << " is empty !";
1062         throw INTERP_KERNEL::Exception(oss.str().c_str());
1063       }
1064   const INTERP_KERNEL::CellModel *cm=&(a[0]->getCellModel());
1065   for(std::size_t ii=0;ii<sz;ii++)
1066     if(&(a[ii]->getCellModel())!=cm)
1067       throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::Merge1SGTUMeshes : all items must have the same geo type !");
1068   std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCoupling1SGTUMesh> > bb(sz);
1069   std::vector< const MEDCoupling1SGTUMesh * > aa(sz);
1070   int spaceDim=-3;
1071   for(std::size_t i=0;i<sz && spaceDim==-3;i++)
1072     {
1073       const MEDCoupling1SGTUMesh *cur=a[i];
1074       const DataArrayDouble *coo=cur->getCoords();
1075       if(coo)
1076         spaceDim=coo->getNumberOfComponents();
1077     }
1078   if(spaceDim==-3)
1079     throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::Merge1SGTUMeshes : no spaceDim specified ! unable to perform merge !");
1080   for(std::size_t i=0;i<sz;i++)
1081     {
1082       bb[i]=a[i]->buildSetInstanceFromThis(spaceDim);
1083       aa[i]=bb[i];
1084     }
1085   return Merge1SGTUMeshesLL(aa);
1086 }
1087
1088 /*!
1089  * \throw If presence of a null instance in the input vector \a a.
1090  * \throw If a is empty
1091  */
1092 MEDCoupling1SGTUMesh *MEDCoupling1SGTUMesh::Merge1SGTUMeshesOnSameCoords(std::vector<const MEDCoupling1SGTUMesh *>& a) throw(INTERP_KERNEL::Exception)
1093 {
1094   if(a.empty())
1095     throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::Merge1SGTUMeshesOnSameCoords : input array must be NON EMPTY !");
1096   std::vector<const MEDCoupling1SGTUMesh *>::const_iterator it=a.begin();
1097   if(!(*it))
1098     throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::Merge1SGTUMeshesOnSameCoords : null instance in the first element of input vector !");
1099   std::vector<const DataArrayInt *> ncs(a.size());
1100   int nbOfCells=(*it)->getNumberOfCells();
1101   const DataArrayDouble *coords=(*it)->getCoords();
1102   const INTERP_KERNEL::CellModel *cm=&((*it)->getCellModel());
1103   int nbNodesPerCell=(*it)->getNumberOfNodesPerCell();
1104   ncs[0]=(*it)->getNodalConnectivity();
1105   it++;
1106   for(int i=1;it!=a.end();i++,it++)
1107     {
1108       if(!(*it))
1109         throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::Merge1SGTUMeshesOnSameCoords : presence of a null instance in the input vector !");
1110       if(cm!=&((*it)->getCellModel()))
1111         throw INTERP_KERNEL::Exception("Geometric types mismatches, Merge1SGTUMeshes impossible !");
1112       (*it)->getNumberOfCells();//to check that all is OK
1113       ncs[i]=(*it)->getNodalConnectivity();
1114       if(coords!=(*it)->getCoords())
1115         throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::Merge1SGTUMeshesOnSameCoords : not lying on same coords !");
1116     }
1117   MEDCouplingAutoRefCountObjectPtr<MEDCoupling1SGTUMesh> ret(new MEDCoupling1SGTUMesh("merge",*cm));
1118   ret->setCoords(coords);
1119   ret->_conn=DataArrayInt::Aggregate(ncs);
1120   return ret.retn();
1121 }
1122
1123 /*!
1124  * Assume that all instances in \a a are non null. If null it leads to a crash. That's why this method is assigned to be low level (LL)
1125  */
1126 MEDCoupling1SGTUMesh *MEDCoupling1SGTUMesh::Merge1SGTUMeshesLL(std::vector<const MEDCoupling1SGTUMesh *>& a) throw(INTERP_KERNEL::Exception)
1127 {
1128   if(a.empty())
1129     throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::Merge1SGTUMeshes : input array must be NON EMPTY !");
1130   std::vector<const MEDCoupling1SGTUMesh *>::const_iterator it=a.begin();
1131   int nbOfCells=(*it)->getNumberOfCells();
1132   const INTERP_KERNEL::CellModel *cm=&((*it)->getCellModel());
1133   int nbNodesPerCell=(*it)->getNumberOfNodesPerCell();
1134   it++;
1135   for(;it!=a.end();it++)
1136     {
1137       if(cm!=&((*it)->getCellModel()))
1138         throw INTERP_KERNEL::Exception("Geometric types mismatches, Merge1SGTUMeshes impossible !");
1139       nbOfCells+=(*it)->getNumberOfCells();
1140     }
1141   std::vector<const MEDCouplingPointSet *> aps(a.size());
1142   std::copy(a.begin(),a.end(),aps.begin());
1143   MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> pts=MergeNodesArray(aps);
1144   MEDCouplingAutoRefCountObjectPtr<MEDCoupling1SGTUMesh> ret(new MEDCoupling1SGTUMesh("merge",*cm));
1145   ret->setCoords(pts);
1146   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> c=DataArrayInt::New();
1147   c->alloc(nbOfCells*nbNodesPerCell,1);
1148   int *cPtr=c->getPointer();
1149   int offset=0;
1150   for(it=a.begin();it!=a.end();it++)
1151     {
1152       int curConnLgth=(*it)->getNodalConnectivityLength();
1153       const int *curC=(*it)->_conn->begin();
1154       cPtr=std::transform(curC,curC+curConnLgth,cPtr,std::bind2nd(std::plus<int>(),offset));
1155       offset+=(*it)->getNumberOfNodes();
1156     }
1157   //
1158   ret->setNodalConnectivity(c);
1159   return ret.retn();
1160 }
1161
1162 MEDCouplingPointSet *MEDCoupling1SGTUMesh::buildPartOfMySelfKeepCoords(const int *begin, const int *end) const
1163 {
1164   int ncell=getNumberOfCells();
1165   MEDCouplingAutoRefCountObjectPtr<MEDCoupling1SGTUMesh> ret(new MEDCoupling1SGTUMesh(getName().c_str(),*_cm));
1166   ret->setCoords(_coords);
1167   std::size_t nbOfElemsRet=std::distance(begin,end);
1168   const int *inConn=_conn->getConstPointer();
1169   int sz=getNumberOfNodesPerCell();
1170   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> connRet=DataArrayInt::New(); connRet->alloc((int)nbOfElemsRet*sz,1);
1171   int *connPtr=connRet->getPointer();
1172   for(const int *work=begin;work!=end;work++,connPtr+=sz)
1173     {
1174       if(*work>=0 && *work<ncell)
1175         std::copy(inConn+(work[0])*sz,inConn+(work[0]+1)*sz,connPtr);
1176       else
1177         {
1178           std::ostringstream oss; oss << "MEDCoupling1SGTUMesh::buildPartOfMySelfKeepCoords : On pos #" << std::distance(begin,work) << " input cell id =" << *work << " should be in [0," << ncell << ") !";
1179           throw INTERP_KERNEL::Exception(oss.str().c_str());
1180         }
1181     }
1182   ret->_conn=connRet;
1183   ret->copyTinyInfoFrom(this);
1184   return ret.retn();
1185 }
1186
1187 MEDCouplingPointSet *MEDCoupling1SGTUMesh::buildPartOfMySelfKeepCoords2(int start, int end, int step) const
1188 {
1189   int ncell=getNumberOfCells();
1190   int nbOfElemsRet=DataArray::GetNumberOfItemGivenBESRelative(start,end,step,"MEDCoupling1SGTUMesh::buildPartOfMySelfKeepCoords2 : ");
1191   MEDCouplingAutoRefCountObjectPtr<MEDCoupling1SGTUMesh> ret(new MEDCoupling1SGTUMesh(getName().c_str(),*_cm));
1192   ret->setCoords(_coords);
1193   const int *inConn=_conn->getConstPointer();
1194   int sz=getNumberOfNodesPerCell();
1195   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> connRet=DataArrayInt::New(); connRet->alloc((int)nbOfElemsRet*sz,1);
1196   int *connPtr=connRet->getPointer();
1197   int curId=start;
1198   for(int i=0;i<nbOfElemsRet;i++,connPtr+=sz,curId+=step)
1199     {
1200       if(curId>=0 && curId<ncell)
1201         std::copy(inConn+curId*sz,inConn+(curId+1)*sz,connPtr);
1202       else
1203         {
1204           std::ostringstream oss; oss << "MEDCoupling1SGTUMesh::buildPartOfMySelfKeepCoords2 : On pos #" << i << " input cell id =" << curId  << " should be in [0," << ncell << ") !";
1205           throw INTERP_KERNEL::Exception(oss.str().c_str());
1206         }
1207     }
1208   ret->_conn=connRet;
1209   ret->copyTinyInfoFrom(this);
1210   return ret.retn();
1211 }
1212
1213 void MEDCoupling1SGTUMesh::computeNodeIdsAlg(std::vector<bool>& nodeIdsInUse) const throw(INTERP_KERNEL::Exception)
1214 {
1215   int sz((int)nodeIdsInUse.size());
1216   int nbCells(getNumberOfCells());
1217   int nbOfNodesPerCell(getNumberOfNodesPerCell());
1218   const int *w(_conn->begin());
1219   for(int i=0;i<nbCells;i++)
1220     for(int j=0;j<nbOfNodesPerCell;j++,w++)
1221       {
1222         if(*w>=0 && *w<sz)
1223           nodeIdsInUse[*w]=true;
1224         else
1225           {
1226             std::ostringstream oss; oss << "MEDCoupling1SGTUMesh::computeNodeIdsAlg : At cell #" << i << " presence of node id #" << *w << " should be in [0," << sz << ") !";
1227             throw INTERP_KERNEL::Exception(oss.str().c_str());
1228           }
1229       }
1230 }
1231
1232 MEDCoupling1SGTUMesh *MEDCoupling1SGTUMesh::buildSetInstanceFromThis(int spaceDim) const throw(INTERP_KERNEL::Exception)
1233 {
1234   MEDCouplingAutoRefCountObjectPtr<MEDCoupling1SGTUMesh> ret(new MEDCoupling1SGTUMesh(getName().c_str(),*_cm));
1235   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> tmp1;
1236   const DataArrayInt *nodalConn(_conn);
1237   if(!nodalConn)
1238     {
1239       tmp1=DataArrayInt::New(); tmp1->alloc(0,1);
1240     }
1241   else
1242     tmp1=_conn;
1243   ret->_conn=tmp1;
1244   if(!_coords)
1245     {
1246       MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> coords=DataArrayDouble::New(); coords->alloc(0,spaceDim);
1247       ret->setCoords(coords);
1248     }
1249   else
1250     ret->setCoords(_coords);
1251   return ret.retn();
1252 }
1253
1254 DataArrayInt *MEDCoupling1SGTUMesh::simplexizePol0() throw(INTERP_KERNEL::Exception)
1255 {
1256   int nbOfCells=getNumberOfCells();
1257   if(getCellModelEnum()!=INTERP_KERNEL::NORM_QUAD4)
1258     return DataArrayInt::Range(0,nbOfCells,1);
1259   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> newConn=DataArrayInt::New(); newConn->alloc(2*3*nbOfCells,1);
1260   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret=DataArrayInt::New(); ret->alloc(2*nbOfCells,1);
1261   const int *c(_conn->begin());
1262   int *retPtr(ret->getPointer()),*newConnPtr(newConn->getPointer());
1263   for(int i=0;i<nbOfCells;i++,c+=4,newConnPtr+=6,retPtr+=2)
1264     {
1265       newConnPtr[0]=c[0]; newConnPtr[1]=c[1]; newConnPtr[2]=c[2];
1266       newConnPtr[3]=c[0]; newConnPtr[4]=c[2]; newConnPtr[5]=c[3];
1267       retPtr[0]=i; retPtr[1]=i;
1268     }
1269   _conn=newConn;
1270   _cm=&INTERP_KERNEL::CellModel::GetCellModel(INTERP_KERNEL::NORM_TRI3);
1271   updateTime();
1272   return ret.retn();
1273 }
1274
1275 DataArrayInt *MEDCoupling1SGTUMesh::simplexizePol1() throw(INTERP_KERNEL::Exception)
1276 {
1277   int nbOfCells=getNumberOfCells();
1278   if(getCellModelEnum()!=INTERP_KERNEL::NORM_QUAD4)
1279     return DataArrayInt::Range(0,nbOfCells,1);
1280   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> newConn=DataArrayInt::New(); newConn->alloc(2*3*nbOfCells,1);
1281   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret=DataArrayInt::New(); ret->alloc(2*nbOfCells,1);
1282   const int *c(_conn->begin());
1283   int *retPtr(ret->getPointer()),*newConnPtr(newConn->getPointer());
1284   for(int i=0;i<nbOfCells;i++,c+=4,newConnPtr+=6,retPtr+=2)
1285     {
1286       newConnPtr[0]=c[0]; newConnPtr[1]=c[1]; newConnPtr[2]=c[3];
1287       newConnPtr[3]=c[1]; newConnPtr[4]=c[2]; newConnPtr[5]=c[3];
1288       retPtr[0]=i; retPtr[1]=i;
1289     }
1290   _conn=newConn;
1291   _cm=&INTERP_KERNEL::CellModel::GetCellModel(INTERP_KERNEL::NORM_TRI3);
1292   updateTime();
1293   return ret.retn();
1294 }
1295
1296 DataArrayInt *MEDCoupling1SGTUMesh::simplexizePlanarFace5() throw(INTERP_KERNEL::Exception)
1297 {
1298   int nbOfCells=getNumberOfCells();
1299   if(getCellModelEnum()!=INTERP_KERNEL::NORM_HEXA8)
1300     return DataArrayInt::Range(0,nbOfCells,1);
1301   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> newConn=DataArrayInt::New(); newConn->alloc(5*4*nbOfCells,1);
1302   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret=DataArrayInt::New(); ret->alloc(5*nbOfCells,1);
1303   const int *c(_conn->begin());
1304   int *retPtr(ret->getPointer()),*newConnPtr(newConn->getPointer());
1305   for(int i=0;i<nbOfCells;i++,c+=8,newConnPtr+=20,retPtr+=5)
1306     {
1307       for(int j=0;j<20;j++)
1308         newConnPtr[j]=c[INTERP_KERNEL::SPLIT_NODES_5_WO[j]];
1309       retPtr[0]=i; retPtr[1]=i; retPtr[2]=i; retPtr[3]=i; retPtr[4]=i;
1310     }
1311   _conn=newConn;
1312   _cm=&INTERP_KERNEL::CellModel::GetCellModel(INTERP_KERNEL::NORM_TETRA4);
1313   updateTime();
1314   return ret.retn();
1315 }
1316
1317 DataArrayInt *MEDCoupling1SGTUMesh::simplexizePlanarFace6() throw(INTERP_KERNEL::Exception)
1318 {
1319   int nbOfCells=getNumberOfCells();
1320   if(getCellModelEnum()!=INTERP_KERNEL::NORM_HEXA8)
1321     return DataArrayInt::Range(0,nbOfCells,1);
1322   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> newConn=DataArrayInt::New(); newConn->alloc(6*4*nbOfCells,1);
1323   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret=DataArrayInt::New(); ret->alloc(6*nbOfCells,1);
1324   const int *c(_conn->begin());
1325   int *retPtr(ret->getPointer()),*newConnPtr(newConn->getPointer());
1326   for(int i=0;i<nbOfCells;i++,c+=8,newConnPtr+=24,retPtr+=6)
1327     {
1328       for(int j=0;j<24;j++)
1329         newConnPtr[j]=c[INTERP_KERNEL::SPLIT_NODES_6_WO[j]];
1330       retPtr[0]=i; retPtr[1]=i; retPtr[2]=i; retPtr[3]=i; retPtr[4]=i; retPtr[5]=i;
1331     }
1332   _conn=newConn;
1333   _cm=&INTERP_KERNEL::CellModel::GetCellModel(INTERP_KERNEL::NORM_TETRA4);
1334   updateTime();
1335   return ret.retn();
1336 }
1337
1338 void MEDCoupling1SGTUMesh::reprQuickOverview(std::ostream& stream) const throw(INTERP_KERNEL::Exception)
1339 {
1340   stream << "MEDCoupling1SGTUMesh C++ instance at " << this << ". Type=" << _cm->getRepr() << ". Name : \"" << getName() << "\".";
1341   stream << " Mesh dimension : " << getMeshDimension() << ".";
1342   if(!_coords)
1343     { stream << " No coordinates set !"; return ; }
1344   if(!_coords->isAllocated())
1345     { stream << " Coordinates set but not allocated !"; return ; }
1346   stream << " Space dimension : " << _coords->getNumberOfComponents() << "." << std::endl;
1347   stream << "Number of nodes : " << _coords->getNumberOfTuples() << ".";
1348   if(!(const DataArrayInt *)_conn)
1349     { stream << std::endl << "Nodal connectivity NOT set !"; return ; }
1350   if(_conn->isAllocated())
1351     {
1352       if(_conn->getNumberOfComponents()==1)
1353         stream << std::endl << "Number of cells : " << getNumberOfCells() << ".";
1354     }
1355 }
1356
1357 void MEDCoupling1SGTUMesh::checkFullyDefined() const throw(INTERP_KERNEL::Exception)
1358 {
1359   if(!((const DataArrayInt *)_conn) || !((const DataArrayDouble *)_coords))
1360     throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::checkFullyDefined : part of this is not fully defined.");
1361 }
1362
1363 /*!
1364  * First step of unserialization process.
1365  */
1366 bool MEDCoupling1SGTUMesh::isEmptyMesh(const std::vector<int>& tinyInfo) const
1367 {
1368   throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::isEmptyMesh : not implemented yet !");
1369 }
1370
1371 void MEDCoupling1SGTUMesh::getTinySerializationInformation(std::vector<double>& tinyInfoD, std::vector<int>& tinyInfo, std::vector<std::string>& littleStrings) const
1372 {
1373   int it,order;
1374   double time=getTime(it,order);
1375   tinyInfo.clear(); tinyInfoD.clear(); littleStrings.clear();
1376   //
1377   littleStrings.push_back(getName());
1378   littleStrings.push_back(getDescription());
1379   littleStrings.push_back(getTimeUnit());
1380   //
1381   std::vector<std::string> littleStrings2,littleStrings3;
1382   if((const DataArrayDouble *)_coords)
1383     _coords->getTinySerializationStrInformation(littleStrings2);
1384   if((const DataArrayInt *)_conn)
1385     _conn->getTinySerializationStrInformation(littleStrings3);
1386   int sz0((int)littleStrings2.size()),sz1((int)littleStrings3.size());
1387   littleStrings.insert(littleStrings.end(),littleStrings2.begin(),littleStrings2.end());
1388   littleStrings.insert(littleStrings.end(),littleStrings3.begin(),littleStrings3.end());
1389   //
1390   tinyInfo.push_back(getCellModelEnum());
1391   tinyInfo.push_back(it);
1392   tinyInfo.push_back(order);
1393   std::vector<int> tinyInfo2,tinyInfo3;
1394   if((const DataArrayDouble *)_coords)
1395     _coords->getTinySerializationIntInformation(tinyInfo2);
1396   if((const DataArrayInt *)_conn)
1397     _conn->getTinySerializationIntInformation(tinyInfo3);
1398   int sz2((int)tinyInfo2.size()),sz3((int)tinyInfo3.size());
1399   tinyInfo.push_back(sz0); tinyInfo.push_back(sz1); tinyInfo.push_back(sz2); tinyInfo.push_back(sz3);
1400   tinyInfo.insert(tinyInfo.end(),tinyInfo2.begin(),tinyInfo2.end());
1401   tinyInfo.insert(tinyInfo.end(),tinyInfo3.begin(),tinyInfo3.end());
1402   //
1403   tinyInfoD.push_back(time);
1404 }
1405
1406 void MEDCoupling1SGTUMesh::resizeForUnserialization(const std::vector<int>& tinyInfo, DataArrayInt *a1, DataArrayDouble *a2, std::vector<std::string>& littleStrings) const
1407 {
1408   std::vector<int> tinyInfo2(tinyInfo.begin()+7,tinyInfo.begin()+7+tinyInfo[5]);
1409   std::vector<int> tinyInfo1(tinyInfo.begin()+7+tinyInfo[5],tinyInfo.begin()+7+tinyInfo[5]+tinyInfo[6]);
1410   a1->resizeForUnserialization(tinyInfo1);
1411   a2->resizeForUnserialization(tinyInfo2);
1412 }
1413
1414 void MEDCoupling1SGTUMesh::serialize(DataArrayInt *&a1, DataArrayDouble *&a2) const
1415 {
1416   int sz(0);
1417   if((const DataArrayInt *)_conn)
1418     if(_conn->isAllocated())
1419       sz=_conn->getNbOfElems();
1420   a1=DataArrayInt::New();
1421   a1->alloc(sz,1);
1422   if(sz!=0 && (const DataArrayInt *)_conn)
1423     std::copy(_conn->begin(),_conn->end(),a1->getPointer());
1424   sz=0;
1425   if((const DataArrayDouble *)_coords)
1426     if(_coords->isAllocated())
1427       sz=_coords->getNbOfElems();
1428   a2=DataArrayDouble::New();
1429   a2->alloc(sz,1);
1430   if(sz!=0 && (const DataArrayDouble *)_coords)
1431     std::copy(_coords->begin(),_coords->end(),a2->getPointer());
1432 }
1433
1434 void MEDCoupling1SGTUMesh::unserialization(const std::vector<double>& tinyInfoD, const std::vector<int>& tinyInfo, const DataArrayInt *a1, DataArrayDouble *a2,
1435                                            const std::vector<std::string>& littleStrings)
1436 {
1437   INTERP_KERNEL::NormalizedCellType gt((INTERP_KERNEL::NormalizedCellType)tinyInfo[0]);
1438   _cm=&INTERP_KERNEL::CellModel::GetCellModel(gt);
1439   setName(littleStrings[0].c_str());
1440   setDescription(littleStrings[1].c_str());
1441   setTimeUnit(littleStrings[2].c_str());
1442   setTime(tinyInfoD[0],tinyInfo[1],tinyInfo[2]);
1443   int sz0(tinyInfo[3]),sz1(tinyInfo[4]),sz2(tinyInfo[5]),sz3(tinyInfo[6]);
1444   //
1445   _coords=DataArrayDouble::New();
1446   std::vector<int> tinyInfo2(tinyInfo.begin()+7,tinyInfo.begin()+7+sz2);
1447   _coords->resizeForUnserialization(tinyInfo2);
1448   std::copy(a2->begin(),a2->end(),_coords->getPointer());
1449   _conn=DataArrayInt::New();
1450   std::vector<int> tinyInfo3(tinyInfo.begin()+7+sz2,tinyInfo.begin()+7+sz2+sz3);
1451   _conn->resizeForUnserialization(tinyInfo3);
1452   std::copy(a1->begin(),a1->end(),_conn->getPointer());
1453   std::vector<std::string> littleStrings2(littleStrings.begin()+3,littleStrings.begin()+3+sz0);
1454   _coords->finishUnserialization(tinyInfo2,littleStrings2);
1455   std::vector<std::string> littleStrings3(littleStrings.begin()+3+sz0,littleStrings.begin()+3+sz0+sz1);
1456   _conn->finishUnserialization(tinyInfo3,littleStrings3);
1457 }
1458
1459 /*!
1460  * Checks if \a this and \a other meshes are geometrically equivalent with high
1461  * probability, else an exception is thrown. The meshes are considered equivalent if
1462  * (1) meshes contain the same number of nodes and the same number of elements of the
1463  * same types (2) three cells of the two meshes (first, last and middle) are based
1464  * on coincident nodes (with a specified precision).
1465  *  \param [in] other - the mesh to compare with.
1466  *  \param [in] prec - the precision used to compare nodes of the two meshes.
1467  *  \throw If the two meshes do not match.
1468  */
1469 void MEDCoupling1SGTUMesh::checkFastEquivalWith(const MEDCouplingMesh *other, double prec) const throw(INTERP_KERNEL::Exception)
1470 {
1471   MEDCouplingPointSet::checkFastEquivalWith(other,prec);
1472   const MEDCoupling1SGTUMesh *otherC=dynamic_cast<const MEDCoupling1SGTUMesh *>(other);
1473   if(!otherC)
1474     throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::checkFastEquivalWith : Two meshes are not unstructured with single static geometric type !");
1475   const DataArrayInt *c1(_conn),*c2(otherC->_conn);
1476   if(c1==c2)
1477     return;
1478   if(!c1 || !c2)
1479     throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::checkFastEquivalWith : presence of nodal connectivity only in one of the 2 meshes !");
1480   if((c1->isAllocated() && !c2->isAllocated()) || (!c1->isAllocated() && c2->isAllocated()))
1481     throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::checkFastEquivalWith : in nodal connectivity, only one is allocated !");
1482   if(c1->getNumberOfComponents()!=1 || c1->getNumberOfComponents()!=1)
1483     throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::checkFastEquivalWith : in nodal connectivity, must have 1 and only 1 component !");
1484   if(c1->getHashCode()!=c2->getHashCode())
1485     throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::checkFastEquivalWith : nodal connectivity differs");
1486 }
1487
1488 MEDCouplingPointSet *MEDCoupling1SGTUMesh::mergeMyselfWithOnSameCoords(const MEDCouplingPointSet *other) const
1489 {
1490   if(!other)
1491     throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::mergeMyselfWithOnSameCoords : input other is null !");
1492   const MEDCoupling1SGTUMesh *otherC=dynamic_cast<const MEDCoupling1SGTUMesh *>(other);
1493   if(!otherC)
1494     throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::mergeMyselfWithOnSameCoords : the input other mesh is not of type single statuc geo type unstructured !");
1495   std::vector<const MEDCoupling1SGTUMesh *> ms(2);
1496   ms[0]=this;
1497   ms[1]=otherC;
1498   return Merge1SGTUMeshesOnSameCoords(ms);
1499 }
1500
1501 void MEDCoupling1SGTUMesh::getReverseNodalConnectivity(DataArrayInt *revNodal, DataArrayInt *revNodalIndx) const throw(INTERP_KERNEL::Exception)
1502 {
1503   checkFullyDefined();
1504   int nbOfNodes=getNumberOfNodes();
1505   int *revNodalIndxPtr=(int *)malloc((nbOfNodes+1)*sizeof(int));
1506   revNodalIndx->useArray(revNodalIndxPtr,true,C_DEALLOC,nbOfNodes+1,1);
1507   std::fill(revNodalIndxPtr,revNodalIndxPtr+nbOfNodes+1,0);
1508   const int *conn=_conn->begin();
1509   int nbOfCells=getNumberOfCells();
1510   int nbOfEltsInRevNodal=0;
1511   int nbOfNodesPerCell=getNumberOfNodesPerCell();
1512   for(int eltId=0;eltId<nbOfCells;eltId++)
1513     {
1514       for(int j=0;j<nbOfNodesPerCell;j++,conn++)
1515         {
1516           if(conn[0]>=0 && conn[0]<nbOfNodes)
1517             {
1518               nbOfEltsInRevNodal++;
1519               revNodalIndxPtr[conn[0]+1]++;
1520             }
1521           else
1522             {
1523               std::ostringstream oss; oss << "MEDCoupling1SGTUMesh::getReverseNodalConnectivity : At cell #" << eltId << " presence of nodeId #" << conn[0] << " should be in [0," << nbOfNodes << ") !";
1524               throw INTERP_KERNEL::Exception(oss.str().c_str());
1525             }
1526         }
1527     }
1528   std::transform(revNodalIndxPtr+1,revNodalIndxPtr+nbOfNodes+1,revNodalIndxPtr,revNodalIndxPtr+1,std::plus<int>());
1529   conn=_conn->begin();
1530   int *revNodalPtr=(int *)malloc((nbOfEltsInRevNodal)*sizeof(int));
1531   revNodal->useArray(revNodalPtr,true,C_DEALLOC,nbOfEltsInRevNodal,1);
1532   std::fill(revNodalPtr,revNodalPtr+nbOfEltsInRevNodal,-1);
1533   for(int eltId=0;eltId<nbOfCells;eltId++)
1534     {
1535       for(int j=0;j<nbOfNodesPerCell;j++,conn++)
1536         {
1537           *std::find_if(revNodalPtr+revNodalIndxPtr[*conn],revNodalPtr+revNodalIndxPtr[*conn+1],std::bind2nd(std::equal_to<int>(),-1))=eltId;
1538         }
1539     }
1540 }
1541
1542 /*!
1543  * Use \a nodalConn array as nodal connectivity of \a this. The input \a nodalConn pointer can be null.
1544  */
1545 void MEDCoupling1SGTUMesh::setNodalConnectivity(DataArrayInt *nodalConn) throw(INTERP_KERNEL::Exception)
1546 {
1547   if(nodalConn)
1548     nodalConn->incrRef();
1549   _conn=nodalConn;
1550   declareAsNew();
1551 }
1552
1553 /*!
1554  * \return DataArrayInt * - the internal reference to the nodal connectivity. The caller is not reponsible to deallocate it.
1555  */
1556 DataArrayInt *MEDCoupling1SGTUMesh::getNodalConnectivity() const throw(INTERP_KERNEL::Exception)
1557 {
1558   const DataArrayInt *ret(_conn);
1559   return const_cast<DataArrayInt *>(ret);
1560 }
1561
1562 /*!
1563  * Allocates memory to store an estimation of the given number of cells. Closer is the estimation to the number of cells effectively inserted,
1564  * less will be the needs to realloc. If the number of cells to be inserted is not known simply put 0 to this parameter.
1565  * If a nodal connectivity previouly existed before the call of this method, it will be reset.
1566  *
1567  *  \param [in] nbOfCells - estimation of the number of cell \a this mesh will contain.
1568  */
1569 void MEDCoupling1SGTUMesh::allocateCells(int nbOfCells) throw(INTERP_KERNEL::Exception)
1570 {
1571   if(nbOfCells<0)
1572     throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::allocateCells : the input number of cells should be >= 0 !");
1573   _conn=DataArrayInt::New();
1574   _conn->reserve(getNumberOfNodesPerCell()*nbOfCells);
1575   declareAsNew();
1576 }
1577
1578 /*!
1579  * Appends at the end of \a this a cell having nodal connectivity array defined in [ \a nodalConnOfCellBg, \a nodalConnOfCellEnd ).
1580  *
1581  * \param [in] nodalConnOfCellBg - the begin (included) of nodal connectivity of the cell to add.
1582  * \param [in] nodalConnOfCellEnd - the end (excluded) of nodal connectivity of the cell to add.
1583  * \throw If the length of the input nodal connectivity array of the cell to add is not equal to number of nodes per cell relative to the unique geometric type
1584  *        attached to \a this.
1585  * \thow If the nodal connectivity array in \a this is null (call MEDCoupling1SGTUMesh::allocateCells before).
1586  */
1587 void MEDCoupling1SGTUMesh::insertNextCell(const int *nodalConnOfCellBg, const int *nodalConnOfCellEnd) throw(INTERP_KERNEL::Exception)
1588 {
1589   int sz=(int)std::distance(nodalConnOfCellBg,nodalConnOfCellEnd);
1590   int ref=getNumberOfNodesPerCell();
1591   if(sz==ref)
1592     {
1593       DataArrayInt *c(_conn);
1594       if(c)
1595         c->pushBackValsSilent(nodalConnOfCellBg,nodalConnOfCellEnd);
1596       else
1597         throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::insertNextCell : nodal connectivity array is null ! Call MEDCoupling1SGTUMesh::allocateCells before !");
1598     }
1599   else
1600     {
1601       std::ostringstream oss; oss << "MEDCoupling1SGTUMesh::insertNextCell : input nodal size (" << sz << ") does not match number of nodes per cell of this (";
1602       oss << ref << ") !";
1603       throw INTERP_KERNEL::Exception(oss.str().c_str());
1604     }
1605 }
1606
1607 /*!
1608  * This method builds the dual mesh of \a this and returns it.
1609  *
1610  * \return MEDCoupling1SGTUMesh * - newly object created to be managed by the caller.
1611  * \throw If \a this is not a mesh containing only simplex cells.
1612  * \throw If \a this is not correctly allocated (coordinates and connectivities have to be correctly set !).
1613  * \throw If at least one node in \a this is orphan (without any simplex cell lying on it !)
1614  */
1615 MEDCoupling1GTUMesh *MEDCoupling1SGTUMesh::computeDualMesh() const throw(INTERP_KERNEL::Exception)
1616 {
1617   const INTERP_KERNEL::CellModel& cm(getCellModel());
1618   if(!cm.isSimplex())
1619     throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::computeDualMesh : this mesh is not a simplex mesh ! Please invoke simplexize of tetrahedrize on this before calling this method !");
1620   switch(getMeshDimension())
1621     {
1622     case 3:
1623       return computeDualMesh3D();
1624     case 2:
1625       return computeDualMesh2D();
1626     default:
1627       throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::computeDualMesh : meshdimension must be in [2,3] !");
1628     }
1629 }
1630
1631 MEDCoupling1DGTUMesh *MEDCoupling1SGTUMesh::computeDualMesh3D() const throw(INTERP_KERNEL::Exception)
1632 {
1633   static const int DUAL_TETRA_0[36]={
1634     4,1,0, 6,0,3, 7,3,1,
1635     4,0,1, 5,2,0, 8,1,2,
1636     6,3,0, 5,0,2, 9,2,3,
1637     7,1,3, 9,3,2, 8,2,1
1638   };
1639   static const int DUAL_TETRA_1[36]={
1640     8,4,10, 11,5,8, 10,7,11,
1641     9,4,8, 8,5,12, 12,6,9,
1642     10,4,9, 9,6,13, 13,7,10,
1643     12,5,11, 13,6,12, 11,7,13
1644   };
1645   static const int FACEID_NOT_SH_NODE[4]={2,3,1,0};
1646   if(getCellModelEnum()!=INTERP_KERNEL::NORM_TETRA4)
1647     throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::computeDualMesh3D : only TETRA4 supported !");
1648   checkFullyDefined();
1649   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> thisu(buildUnstructured());
1650   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> revNodArr(DataArrayInt::New()),revNodIArr(DataArrayInt::New());
1651   thisu->getReverseNodalConnectivity(revNodArr,revNodIArr);
1652   const int *revNod(revNodArr->begin()),*revNodI(revNodIArr->begin()),*nodal(_conn->begin());
1653   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> d1Arr(DataArrayInt::New()),di1Arr(DataArrayInt::New()),rd1Arr(DataArrayInt::New()),rdi1Arr(DataArrayInt::New());
1654   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> edges(thisu->explode3DMeshTo1D(d1Arr,di1Arr,rd1Arr,rdi1Arr));
1655   const int *d1(d1Arr->begin());
1656   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> d2Arr(DataArrayInt::New()),di2Arr(DataArrayInt::New()),rd2Arr(DataArrayInt::New()),rdi2Arr(DataArrayInt::New());
1657   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> faces(thisu->buildDescendingConnectivity(d2Arr,di2Arr,rd2Arr,rdi2Arr));  thisu=0;
1658   const int *d2(d2Arr->begin()),*rd2(rd2Arr->begin()),*rdi2(rdi2Arr->begin());
1659   MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> edgesBaryArr(edges->getBarycenterAndOwner()),facesBaryArr(faces->getBarycenterAndOwner()),baryArr(getBarycenterAndOwner());
1660   const int nbOfNodes(getNumberOfNodes()),offset0(nbOfNodes+faces->getNumberOfCells()),offset1(offset0+edges->getNumberOfCells());
1661   edges=0; faces=0;
1662   std::vector<const DataArrayDouble *> v(4); v[0]=getCoords(); v[1]=facesBaryArr; v[2]=edgesBaryArr; v[3]=baryArr;
1663   MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> zeArr(DataArrayDouble::Aggregate(v)); baryArr=0; edgesBaryArr=0; facesBaryArr=0;
1664   std::string name("DualOf_"); name+=getName();
1665   MEDCouplingAutoRefCountObjectPtr<MEDCoupling1DGTUMesh> ret(MEDCoupling1DGTUMesh::New(name.c_str(),INTERP_KERNEL::NORM_POLYHED)); ret->setCoords(zeArr);
1666   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> cArr(DataArrayInt::New()),ciArr(DataArrayInt::New()); ciArr->alloc(nbOfNodes+1,1); ciArr->setIJ(0,0,0); cArr->alloc(0,1);
1667   for(int i=0;i<nbOfNodes;i++,revNodI++)
1668     {
1669       int nbOfCellsSharingNode(revNodI[1]-revNodI[0]);
1670       if(nbOfCellsSharingNode==0)
1671         {
1672           std::ostringstream oss; oss << "MEDCoupling1SGTUMesh::computeDualMesh3D : Node #" << i << " is orphan !";
1673           throw INTERP_KERNEL::Exception(oss.str().c_str());
1674         }
1675       for(int j=0;j<nbOfCellsSharingNode;j++)
1676         {
1677           int curCellId(revNod[revNodI[0]+j]);
1678           const int *connOfCurCell(nodal+4*curCellId);
1679           std::size_t nodePosInCurCell(std::distance(connOfCurCell,std::find(connOfCurCell,connOfCurCell+4,i)));
1680           if(j!=0) cArr->pushBackSilent(-1);
1681           int tmp[14];
1682           //
1683           tmp[0]=d1[6*curCellId+DUAL_TETRA_0[nodePosInCurCell*9+0]-4]+offset0; tmp[1]=d2[4*curCellId+DUAL_TETRA_0[nodePosInCurCell*9+1]]+nbOfNodes;
1684           tmp[2]=curCellId+offset1; tmp[3]=d2[4*curCellId+DUAL_TETRA_0[nodePosInCurCell*9+2]]+nbOfNodes;
1685           tmp[4]=-1;
1686           tmp[5]=d1[6*curCellId+DUAL_TETRA_0[nodePosInCurCell*9+3]-4]+offset0; tmp[6]=d2[4*curCellId+DUAL_TETRA_0[nodePosInCurCell*9+4]]+nbOfNodes;
1687           tmp[7]=curCellId+offset1; tmp[8]=d2[4*curCellId+DUAL_TETRA_0[nodePosInCurCell*9+5]]+nbOfNodes;
1688           tmp[9]=-1;
1689           tmp[10]=d1[6*curCellId+DUAL_TETRA_0[nodePosInCurCell*9+6]-4]+offset0; tmp[11]=d2[4*curCellId+DUAL_TETRA_0[nodePosInCurCell*9+7]]+nbOfNodes;
1690           tmp[12]=curCellId+offset1; tmp[13]=d2[4*curCellId+DUAL_TETRA_0[nodePosInCurCell*9+8]]+nbOfNodes;
1691           cArr->insertAtTheEnd(tmp,tmp+14);
1692           int kk(0);
1693           for(int k=0;k<4;k++)
1694             {
1695               if(FACEID_NOT_SH_NODE[nodePosInCurCell]!=k)
1696                 {
1697                   const int *faceId(d2+4*curCellId+k);
1698                   if(rdi2[*faceId+1]-rdi2[*faceId]==1)
1699                     {
1700                       int tmp2[5]; tmp2[0]=-1; tmp2[1]=i;
1701                       tmp2[2]=d1[6*curCellId+DUAL_TETRA_1[9*nodePosInCurCell+3*kk+0]-8]+offset0;
1702                       tmp2[3]=d2[4*curCellId+DUAL_TETRA_1[9*nodePosInCurCell+3*kk+1]-4]+nbOfNodes;
1703                       tmp2[4]=d1[6*curCellId+DUAL_TETRA_1[9*nodePosInCurCell+3*kk+2]-8]+offset0;
1704                       cArr->insertAtTheEnd(tmp2,tmp2+5);
1705                     }
1706                   kk++;
1707                 }
1708             }
1709         }
1710       ciArr->setIJ(i+1,0,cArr->getNumberOfTuples());
1711     }
1712   ret->setNodalConnectivity(cArr,ciArr);
1713   return ret.retn();
1714 }
1715
1716 MEDCoupling1DGTUMesh *MEDCoupling1SGTUMesh::computeDualMesh2D() const throw(INTERP_KERNEL::Exception)
1717 {
1718   static const int DUAL_TRI_0[6]={0,2, 1,0, 2,1};
1719   static const int DUAL_TRI_1[6]={-3,+5, +3,-4, +4,-5};
1720   static const int FACEID_NOT_SH_NODE[3]={1,2,0};
1721   if(getCellModelEnum()!=INTERP_KERNEL::NORM_TRI3)
1722     throw INTERP_KERNEL::Exception("MEDCoupling1SGTUMesh::computeDualMesh2D : only TRI3 supported !");
1723   checkFullyDefined();
1724   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> thisu(buildUnstructured());
1725   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> revNodArr(DataArrayInt::New()),revNodIArr(DataArrayInt::New());
1726   thisu->getReverseNodalConnectivity(revNodArr,revNodIArr);
1727   const int *revNod(revNodArr->begin()),*revNodI(revNodIArr->begin()),*nodal(_conn->begin());
1728   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> d2Arr(DataArrayInt::New()),di2Arr(DataArrayInt::New()),rd2Arr(DataArrayInt::New()),rdi2Arr(DataArrayInt::New());
1729   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> edges(thisu->buildDescendingConnectivity(d2Arr,di2Arr,rd2Arr,rdi2Arr));  thisu=0;
1730   const int *d2(d2Arr->begin()),*rd2(rd2Arr->begin()),*rdi2(rdi2Arr->begin());
1731   MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> edgesBaryArr(edges->getBarycenterAndOwner()),baryArr(getBarycenterAndOwner());
1732   const int nbOfNodes(getNumberOfNodes()),offset0(nbOfNodes+edges->getNumberOfCells());
1733   edges=0;
1734   std::vector<const DataArrayDouble *> v(3); v[0]=getCoords(); v[1]=edgesBaryArr; v[2]=baryArr;
1735   MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> zeArr(DataArrayDouble::Aggregate(v)); baryArr=0; edgesBaryArr=0;
1736   std::string name("DualOf_"); name+=getName();
1737   MEDCouplingAutoRefCountObjectPtr<MEDCoupling1DGTUMesh> ret(MEDCoupling1DGTUMesh::New(name.c_str(),INTERP_KERNEL::NORM_POLYGON)); ret->setCoords(zeArr);
1738   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> cArr(DataArrayInt::New()),ciArr(DataArrayInt::New()); ciArr->alloc(nbOfNodes+1,1); ciArr->setIJ(0,0,0); cArr->alloc(0,1);
1739   for(int i=0;i<nbOfNodes;i++,revNodI++)
1740     {
1741       int nbOfCellsSharingNode(revNodI[1]-revNodI[0]);
1742       if(nbOfCellsSharingNode==0)
1743         {
1744           std::ostringstream oss; oss << "MEDCoupling1SGTUMesh::computeDualMesh2D : Node #" << i << " is orphan !";
1745           throw INTERP_KERNEL::Exception(oss.str().c_str());
1746         }
1747       std::vector< std::vector<int> > polyg;
1748       for(int j=0;j<nbOfCellsSharingNode;j++)
1749         {
1750           int curCellId(revNod[revNodI[0]+j]);
1751           const int *connOfCurCell(nodal+3*curCellId);
1752           std::size_t nodePosInCurCell(std::distance(connOfCurCell,std::find(connOfCurCell,connOfCurCell+4,i)));
1753           std::vector<int> locV(3);
1754           locV[0]=d2[3*curCellId+DUAL_TRI_0[2*nodePosInCurCell+0]]+nbOfNodes; locV[1]=curCellId+offset0; locV[2]=d2[3*curCellId+DUAL_TRI_0[2*nodePosInCurCell+1]]+nbOfNodes;
1755           polyg.push_back(locV);
1756           int kk(0);
1757           for(int k=0;k<3;k++)
1758             {
1759               if(FACEID_NOT_SH_NODE[nodePosInCurCell]!=k)
1760                 {
1761                   const int *edgeId(d2+3*curCellId+k);
1762                   if(rdi2[*edgeId+1]-rdi2[*edgeId]==1)
1763                     {
1764                       std::vector<int> locV2(2);
1765                       int zeLocEdgeIdRel(DUAL_TRI_1[2*nodePosInCurCell+kk]);
1766                       if(zeLocEdgeIdRel>0)
1767                         {  locV2[0]=d2[3*curCellId+zeLocEdgeIdRel-3]+nbOfNodes;  locV2[1]=i; }
1768                       else
1769                         {  locV2[0]=i; locV2[1]=d2[3*curCellId-zeLocEdgeIdRel-3]+nbOfNodes; }
1770                       polyg.push_back(locV2);
1771                     }
1772                   kk++;
1773                 }
1774             }
1775         }
1776       std::vector<int> zePolyg(MEDCoupling1DGTUMesh::BuildAPolygonFromParts(polyg));
1777       cArr->insertAtTheEnd(zePolyg.begin(),zePolyg.end());
1778       ciArr->setIJ(i+1,0,cArr->getNumberOfTuples());
1779     }
1780   ret->setNodalConnectivity(cArr,ciArr);
1781   return ret.retn();
1782 }
1783
1784 /*!
1785  * This method aggregate the bbox of each cell and put it into bbox
1786  *
1787  * \return DataArrayDouble * - newly created object (to be managed by the caller) \a this number of cells tuples and 2*spacedim components.
1788  *
1789  * \throw If \a this is not fully set (coordinates and connectivity).
1790  * \throw If a cell in \a this has no valid nodeId.
1791  */
1792 DataArrayDouble *MEDCoupling1SGTUMesh::getBoundingBoxForBBTree() const
1793 {
1794   int spaceDim(getSpaceDimension()),nbOfCells(getNumberOfCells()),nbOfNodes(getNumberOfNodes()),nbOfNodesPerCell(getNumberOfNodesPerCell());
1795   MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> ret(DataArrayDouble::New()); ret->alloc(nbOfCells,2*spaceDim);
1796   double *bbox(ret->getPointer());
1797   for(int i=0;i<nbOfCells*spaceDim;i++)
1798     {
1799       bbox[2*i]=std::numeric_limits<double>::max();
1800       bbox[2*i+1]=-std::numeric_limits<double>::max();
1801     }
1802   const double *coordsPtr(_coords->getConstPointer());
1803   const int *conn(_conn->getConstPointer());
1804   for(int i=0;i<nbOfCells;i++)
1805     {
1806       int kk(0);
1807       for(int j=0;j<nbOfNodesPerCell;j++,conn++)
1808         {
1809           int nodeId(*conn);
1810           if(nodeId>=0 && nodeId<nbOfNodes)
1811             {
1812               for(int k=0;k<spaceDim;k++)
1813                 {
1814                   bbox[2*spaceDim*i+2*k]=std::min(bbox[2*spaceDim*i+2*k],coordsPtr[spaceDim*nodeId+k]);
1815                   bbox[2*spaceDim*i+2*k+1]=std::max(bbox[2*spaceDim*i+2*k+1],coordsPtr[spaceDim*nodeId+k]);
1816                 }
1817               kk++;
1818             }
1819         }
1820       if(kk==0)
1821         {
1822           std::ostringstream oss; oss << "MEDCoupling1SGTUMesh::getBoundingBoxForBBTree : cell #" << i << " contains no valid nodeId !";
1823           throw INTERP_KERNEL::Exception(oss.str().c_str());
1824         }
1825     }
1826   return ret.retn();
1827 }
1828
1829 //==
1830
1831 MEDCoupling1DGTUMesh *MEDCoupling1DGTUMesh::New()
1832 {
1833   return new MEDCoupling1DGTUMesh;
1834 }
1835
1836 MEDCoupling1DGTUMesh *MEDCoupling1DGTUMesh::New(const char *name, INTERP_KERNEL::NormalizedCellType type) throw(INTERP_KERNEL::Exception)
1837 {
1838   if(type==INTERP_KERNEL::NORM_ERROR)
1839     throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::New : NORM_ERROR is not a valid type to be used as base geometric type for a mesh !");
1840   const INTERP_KERNEL::CellModel& cm=INTERP_KERNEL::CellModel::GetCellModel(type);
1841   if(!cm.isDynamic())
1842     {
1843       std::ostringstream oss; oss << "MEDCoupling1DGTUMesh::New : the input geometric type " << cm.getRepr() << " is static ! Only dynamic types are allowed here !";
1844       throw INTERP_KERNEL::Exception(oss.str().c_str());
1845     }
1846   return new MEDCoupling1DGTUMesh(name,cm);
1847 }
1848
1849 MEDCoupling1DGTUMesh::MEDCoupling1DGTUMesh()
1850 {
1851 }
1852
1853 MEDCoupling1DGTUMesh::MEDCoupling1DGTUMesh(const char *name, const INTERP_KERNEL::CellModel& cm):MEDCoupling1GTUMesh(name,cm)
1854 {
1855 }
1856
1857 MEDCoupling1DGTUMesh::MEDCoupling1DGTUMesh(const MEDCoupling1DGTUMesh& other, bool recDeepCpy):MEDCoupling1GTUMesh(other,recDeepCpy),_conn(other._conn)
1858 {
1859   if(recDeepCpy)
1860     {
1861       const DataArrayInt *c(other._conn);
1862       if(c)
1863         _conn=c->deepCpy();
1864       c=other._conn_indx;
1865       if(c)
1866         _conn_indx=c->deepCpy();
1867     }
1868 }
1869
1870 MEDCoupling1DGTUMesh *MEDCoupling1DGTUMesh::clone(bool recDeepCpy) const
1871 {
1872   return new MEDCoupling1DGTUMesh(*this,recDeepCpy);
1873 }
1874
1875 /*!
1876  * This method behaves mostly like MEDCoupling1DGTUMesh::deepCpy method, except that only nodal connectivity arrays are deeply copied.
1877  * The coordinates are shared between \a this and the returned instance.
1878  *
1879  * \return MEDCouplingUMesh * - A new object instance holding the copy of \a this (deep for connectivity, shallow for coordiantes)
1880  * \sa MEDCoupling1DGTUMesh::deepCpy
1881  */
1882 MEDCouplingPointSet *MEDCoupling1DGTUMesh::deepCpyConnectivityOnly() const throw(INTERP_KERNEL::Exception)
1883 {
1884   checkCoherency();
1885   MEDCouplingAutoRefCountObjectPtr<MEDCoupling1DGTUMesh> ret(clone(false));
1886   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> c(_conn->deepCpy()),ci(_conn_indx->deepCpy());
1887   ret->setNodalConnectivity(c,ci);
1888   return ret.retn();
1889 }
1890
1891 void MEDCoupling1DGTUMesh::updateTime() const
1892 {
1893   MEDCoupling1GTUMesh::updateTime();
1894   const DataArrayInt *c(_conn);
1895   if(c)
1896     updateTimeWith(*c);
1897   c=_conn_indx;
1898   if(c)
1899     updateTimeWith(*c);
1900 }
1901
1902 std::size_t MEDCoupling1DGTUMesh::getHeapMemorySize() const
1903 {
1904   std::size_t ret=0;
1905   const DataArrayInt *c(_conn);
1906   if(c)
1907     ret+=c->getHeapMemorySize();
1908   c=_conn_indx;
1909   if(c)
1910     ret+=c->getHeapMemorySize();
1911   return MEDCoupling1GTUMesh::getHeapMemorySize()+ret;
1912 }
1913
1914 MEDCouplingMesh *MEDCoupling1DGTUMesh::deepCpy() const
1915 {
1916   return clone(true);
1917 }
1918
1919 bool MEDCoupling1DGTUMesh::isEqualIfNotWhy(const MEDCouplingMesh *other, double prec, std::string& reason) const throw(INTERP_KERNEL::Exception)
1920 {
1921   if(!other)
1922     throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::isEqualIfNotWhy : input other pointer is null !");
1923   std::ostringstream oss; oss.precision(15);
1924   const MEDCoupling1DGTUMesh *otherC=dynamic_cast<const MEDCoupling1DGTUMesh *>(other);
1925   if(!otherC)
1926     {
1927       reason="mesh given in input is not castable in MEDCoupling1DGTUMesh !";
1928       return false;
1929     }
1930   if(!MEDCoupling1GTUMesh::isEqualIfNotWhy(other,prec,reason))
1931     return false;
1932   const DataArrayInt *c1(_conn),*c2(otherC->_conn);
1933   if(c1==c2)
1934     return true;
1935   if(!c1 || !c2)
1936     {
1937       reason="in connectivity of single dynamic geometric type exactly one among this and other is null !";
1938       return false;
1939     }
1940   if(!c1->isEqualIfNotWhy(*c2,reason))
1941     {
1942       reason.insert(0,"Nodal connectivity DataArrayInt differs : ");
1943       return false;
1944     }
1945   c1=_conn_indx; c2=otherC->_conn_indx;
1946   if(c1==c2)
1947     return true;
1948   if(!c1 || !c2)
1949     {
1950       reason="in connectivity index of single dynamic geometric type exactly one among this and other is null !";
1951       return false;
1952     }
1953   if(!c1->isEqualIfNotWhy(*c2,reason))
1954     {
1955       reason.insert(0,"Nodal connectivity index DataArrayInt differs : ");
1956       return false;
1957     }
1958   return true;
1959 }
1960
1961 bool MEDCoupling1DGTUMesh::isEqualWithoutConsideringStr(const MEDCouplingMesh *other, double prec) const
1962 {
1963   if(!other)
1964     throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::isEqualWithoutConsideringStr : input other pointer is null !");
1965   const MEDCoupling1DGTUMesh *otherC=dynamic_cast<const MEDCoupling1DGTUMesh *>(other);
1966   if(!otherC)
1967     return false;
1968   if(!MEDCoupling1GTUMesh::isEqualWithoutConsideringStr(other,prec))
1969     return false;
1970   const DataArrayInt *c1(_conn),*c2(otherC->_conn);
1971   if(c1==c2)
1972     return true;
1973   if(!c1 || !c2)
1974     return false;
1975   if(!c1->isEqualWithoutConsideringStr(*c2))
1976     return false;
1977   return true;
1978   c1=_conn_indx; c2=otherC->_conn_indx;
1979   if(c1==c2)
1980     return true;
1981   if(!c1 || !c2)
1982     return false;
1983   if(!c1->isEqualWithoutConsideringStr(*c2))
1984     return false;
1985   return true;
1986 }
1987
1988 /*!
1989  * Checks if \a this and \a other meshes are geometrically equivalent with high
1990  * probability, else an exception is thrown. The meshes are considered equivalent if
1991  * (1) meshes contain the same number of nodes and the same number of elements of the
1992  * same types (2) three cells of the two meshes (first, last and middle) are based
1993  * on coincident nodes (with a specified precision).
1994  *  \param [in] other - the mesh to compare with.
1995  *  \param [in] prec - the precision used to compare nodes of the two meshes.
1996  *  \throw If the two meshes do not match.
1997  */
1998 void MEDCoupling1DGTUMesh::checkFastEquivalWith(const MEDCouplingMesh *other, double prec) const throw(INTERP_KERNEL::Exception)
1999 {
2000   MEDCouplingPointSet::checkFastEquivalWith(other,prec);
2001   const MEDCoupling1DGTUMesh *otherC=dynamic_cast<const MEDCoupling1DGTUMesh *>(other);
2002   if(!otherC)
2003     throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::checkFastEquivalWith : Two meshes are not unstructured with single dynamic geometric type !");
2004   const DataArrayInt *c1(_conn),*c2(otherC->_conn);
2005   if(c1!=c2)
2006     {
2007       if(!c1 || !c2)
2008         throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::checkFastEquivalWith : presence of nodal connectivity only in one of the 2 meshes !");
2009       if((c1->isAllocated() && !c2->isAllocated()) || (!c1->isAllocated() && c2->isAllocated()))
2010         throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::checkFastEquivalWith : in nodal connectivity, only one is allocated !");
2011       if(c1->getNumberOfComponents()!=1 || c1->getNumberOfComponents()!=1)
2012         throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::checkFastEquivalWith : in nodal connectivity, must have 1 and only 1 component !");
2013       if(c1->getHashCode()!=c2->getHashCode())
2014         throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::checkFastEquivalWith : nodal connectivity differs");
2015     }
2016   c1=_conn_indx; c2=otherC->_conn_indx;
2017   if(c1!=c2)
2018     {
2019       if(!c1 || !c2)
2020         throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::checkFastEquivalWith : presence of nodal connectivity index only in one of the 2 meshes !");
2021       if((c1->isAllocated() && !c2->isAllocated()) || (!c1->isAllocated() && c2->isAllocated()))
2022         throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::checkFastEquivalWith : in nodal connectivity index, only one is allocated !");
2023       if(c1->getNumberOfComponents()!=1 || c1->getNumberOfComponents()!=1)
2024         throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::checkFastEquivalWith : in nodal connectivity index, must have 1 and only 1 component !");
2025       if(c1->getHashCode()!=c2->getHashCode())
2026         throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::checkFastEquivalWith : nodal connectivity index differs");
2027     }
2028 }
2029
2030 void MEDCoupling1DGTUMesh::checkCoherencyOfConnectivity() const throw(INTERP_KERNEL::Exception)
2031 {
2032   const DataArrayInt *c1(_conn);
2033   if(c1)
2034     {
2035       if(c1->getNumberOfComponents()!=1)
2036         throw INTERP_KERNEL::Exception("Nodal connectivity array is expected to be with number of components set to one !");
2037       if(c1->getInfoOnComponent(0)!="")
2038         throw INTERP_KERNEL::Exception("Nodal connectivity array is expected to have no info on its single component !");
2039       c1->checkAllocated();
2040     }
2041   else
2042     throw INTERP_KERNEL::Exception("Nodal connectivity array not defined !");
2043   //
2044   int sz2=_conn->getNumberOfTuples();
2045   c1=_conn_indx;
2046   if(c1)
2047     {
2048       if(c1->getNumberOfComponents()!=1)
2049         throw INTERP_KERNEL::Exception("Nodal connectivity index array is expected to be with number of components set to one !");
2050       c1->checkAllocated();
2051       if(c1->getNumberOfTuples()<1)
2052         throw INTERP_KERNEL::Exception("Nodal connectivity index array is expected to have a a size of 1 at least !");
2053       if(c1->getInfoOnComponent(0)!="")
2054         throw INTERP_KERNEL::Exception("Nodal connectivity index array is expected to have no info on its single component !");
2055       int f=c1->front(),ll=c1->back();
2056       if(f<0 || f>=sz2)
2057         {
2058           std::ostringstream oss; oss << "Nodal connectivity index array first value (" << f << ") is expected to be exactly in [0," << sz2 << ") !";
2059           throw INTERP_KERNEL::Exception(oss.str().c_str());
2060         }
2061       if(ll<0 || ll>sz2)
2062         {
2063           std::ostringstream oss; oss << "Nodal connectivity index array last value (" << ll << ") is expected to be exactly in [0," << sz2 << "] !";
2064           throw INTERP_KERNEL::Exception(oss.str().c_str());
2065         }
2066       if(f>ll)
2067         {
2068           std::ostringstream oss; oss << "Nodal connectivity index array looks very bad (not increasing monotonic) because front (" << f << ") is greater that back (" << ll << ") !";
2069           throw INTERP_KERNEL::Exception(oss.str().c_str());
2070         }
2071     }
2072   else
2073     throw INTERP_KERNEL::Exception("Nodal connectivity index array not defined !");
2074   int szOfC1Exp=_conn_indx->back();
2075   if(sz2<szOfC1Exp)
2076     {
2077       std::ostringstream oss; oss << "MEDCoupling1DGTUMesh::checkCoherencyOfConnectivity : The expected length of nodal connectivity array regarding index is " << szOfC1Exp << " but the actual size of it is " << c1->getNumberOfTuples() << " !";
2078       throw INTERP_KERNEL::Exception(oss.str().c_str());
2079     }
2080 }
2081
2082 /*!
2083  * If \a this pass this method, you are sure that connectivity arrays are not null, with exactly one component, no name, no component name, allocated.
2084  * In addition you are sure that the length of nodal connectivity index array is bigger than or equal to one.
2085  * In addition you are also sure that length of nodal connectivity is coherent with the content of the last value in the index array.
2086  */
2087 void MEDCoupling1DGTUMesh::checkCoherency() const throw(INTERP_KERNEL::Exception)
2088 {
2089   MEDCouplingPointSet::checkCoherency();
2090   checkCoherencyOfConnectivity();
2091 }
2092
2093 void MEDCoupling1DGTUMesh::checkCoherency1(double eps) const throw(INTERP_KERNEL::Exception)
2094 {
2095   checkCoherency();
2096   const DataArrayInt *c1(_conn),*c2(_conn_indx);
2097   if(!c2->isMonotonic(true))
2098     throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::checkCoherency1 : the nodal connectivity index is expected to be increasing monotinic !");
2099   //
2100   int nbOfTuples=c1->getNumberOfTuples();
2101   int nbOfNodes=getNumberOfNodes();
2102   const int *w(c1->begin());
2103   for(int i=0;i<nbOfTuples;i++,w++)
2104     {
2105       if(*w==-1) continue;
2106       if(*w<0 || *w>=nbOfNodes)
2107         {
2108           std::ostringstream oss; oss << "At pos #" << i << " of nodal connectivity array references to node id #" << *w << " must be in [0," << nbOfNodes << ") !";
2109           throw INTERP_KERNEL::Exception(oss.str().c_str());
2110         }
2111     }
2112 }
2113
2114 void MEDCoupling1DGTUMesh::checkCoherency2(double eps) const throw(INTERP_KERNEL::Exception)
2115 {
2116   checkCoherency1(eps);
2117 }
2118
2119 int MEDCoupling1DGTUMesh::getNumberOfCells() const
2120 {
2121   checkCoherencyOfConnectivity();//do not remove
2122   return _conn_indx->getNumberOfTuples()-1;
2123 }
2124
2125 /*!
2126  * This method returns a newly allocated array containing this->getNumberOfCells() tuples and 1 component.
2127  * For each cell in \b this the number of nodes constituting cell is computed.
2128  * For each polyhedron cell, the sum of the number of nodes of each face constituting polyhedron cell is returned.
2129  * So for pohyhedrons some nodes can be counted several times in the returned result.
2130  *
2131  * \return a newly allocated array
2132  */
2133 DataArrayInt *MEDCoupling1DGTUMesh::computeNbOfNodesPerCell() const throw(INTERP_KERNEL::Exception)
2134 {
2135   checkCoherency();
2136   _conn_indx->checkMonotonic(true);
2137   if(getCellModelEnum()!=INTERP_KERNEL::NORM_POLYHED)
2138     return _conn_indx->deltaShiftIndex();
2139   // for polyhedrons
2140   int nbOfCells=_conn_indx->getNumberOfTuples()-1;
2141   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret=DataArrayInt::New();
2142   ret->alloc(nbOfCells,1);
2143   int *retPtr=ret->getPointer();
2144   const int *ci=_conn_indx->begin(),*c=_conn->begin();
2145   for(int i=0;i<nbOfCells;i++,retPtr++,ci++)
2146     *retPtr=ci[1]-ci[0]-std::count(c+ci[0],c+ci[1],-1);
2147   return ret.retn();
2148 }
2149
2150 /*!
2151  * This method returns a newly allocated array containing this->getNumberOfCells() tuples and 1 component.
2152  * For each cell in \b this the number of faces constituting (entity of dimension this->getMeshDimension()-1) cell is computed.
2153  *
2154  * \return a newly allocated array
2155  */
2156 DataArrayInt *MEDCoupling1DGTUMesh::computeNbOfFacesPerCell() const throw(INTERP_KERNEL::Exception)
2157 {
2158   checkCoherency();
2159   _conn_indx->checkMonotonic(true);
2160   if(getCellModelEnum()!=INTERP_KERNEL::NORM_POLYHED && getCellModelEnum()!=INTERP_KERNEL::NORM_QPOLYG)
2161     return _conn_indx->deltaShiftIndex();
2162   if(getCellModelEnum()==INTERP_KERNEL::NORM_QPOLYG)
2163     {
2164       MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret=_conn_indx->deltaShiftIndex();
2165       ret->applyDivideBy(2);
2166       return ret.retn();
2167     }
2168   // for polyhedrons
2169   int nbOfCells=_conn_indx->getNumberOfTuples()-1;
2170   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret=DataArrayInt::New();
2171   ret->alloc(nbOfCells,1);
2172   int *retPtr=ret->getPointer();
2173   const int *ci=_conn_indx->begin(),*c=_conn->begin();
2174   for(int i=0;i<nbOfCells;i++,retPtr++,ci++)
2175     *retPtr=std::count(c+ci[0],c+ci[1],-1)+1;
2176   return ret.retn();
2177 }
2178
2179 /*!
2180  * This method computes effective number of nodes per cell. That is to say nodes appearing several times in nodal connectivity of a cell,
2181  * will be counted only once here whereas it will be counted several times in MEDCoupling1DGTUMesh::computeNbOfNodesPerCell method.
2182  *
2183  * \return DataArrayInt * - new object to be deallocated by the caller.
2184  * \sa MEDCoupling1DGTUMesh::computeNbOfNodesPerCell
2185  */
2186 DataArrayInt *MEDCoupling1DGTUMesh::computeEffectiveNbOfNodesPerCell() const throw(INTERP_KERNEL::Exception)
2187 {
2188   checkCoherency();
2189   _conn_indx->checkMonotonic(true);
2190   int nbOfCells(_conn_indx->getNumberOfTuples()-1);
2191   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret=DataArrayInt::New();
2192   ret->alloc(nbOfCells,1);
2193   int *retPtr(ret->getPointer());
2194   const int *ci(_conn_indx->begin()),*c(_conn->begin());
2195   if(getCellModelEnum()!=INTERP_KERNEL::NORM_POLYHED)
2196     {
2197       for(int i=0;i<nbOfCells;i++,retPtr++,ci++)
2198         {
2199           std::set<int> s(c+ci[0],c+ci[1]);
2200           *retPtr=(int)s.size();
2201         }
2202     }
2203   else
2204     {
2205       for(int i=0;i<nbOfCells;i++,retPtr++,ci++)
2206         {
2207           std::set<int> s(c+ci[0],c+ci[1]); s.erase(-1);
2208           *retPtr=(int)s.size();
2209         }
2210     }
2211   return ret.retn();
2212 }
2213
2214 void MEDCoupling1DGTUMesh::getNodeIdsOfCell(int cellId, std::vector<int>& conn) const
2215 {
2216   int nbOfCells(getNumberOfCells());//performs checks
2217   if(cellId>=0 && cellId<nbOfCells)
2218     {
2219       int strt=_conn_indx->getIJ(cellId,0),stp=_conn_indx->getIJ(cellId+1,0);
2220       int nbOfNodes=stp-strt;
2221       if(nbOfNodes<0)
2222         throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::getNodeIdsOfCell : the index array is invalid ! Should be increasing monotonic !");
2223       conn.resize(nbOfNodes);
2224       std::copy(_conn->begin()+strt,_conn->begin()+stp,conn.begin());
2225     }
2226   else
2227     {
2228       std::ostringstream oss; oss << "MEDCoupling1DGTUMesh::getNodeIdsOfCell : request for cellId #" << cellId << " must be in [0," << nbOfCells << ") !";
2229       throw INTERP_KERNEL::Exception(oss.str().c_str());
2230     }
2231 }
2232
2233 int MEDCoupling1DGTUMesh::getNumberOfNodesInCell(int cellId) const throw(INTERP_KERNEL::Exception)
2234 {
2235   int nbOfCells(getNumberOfCells());//performs checks
2236   if(cellId>=0 && cellId<nbOfCells)
2237     {
2238       const int *conn(_conn->begin());
2239       int strt=_conn_indx->getIJ(cellId,0),stp=_conn_indx->getIJ(cellId+1,0);
2240       return stp-strt-std::count(conn+strt,conn+stp,-1);
2241     }
2242   else
2243     {
2244       std::ostringstream oss; oss << "MEDCoupling1DGTUMesh::getNumberOfNodesInCell : request for cellId #" << cellId << " must be in [0," << nbOfCells << ") !";
2245       throw INTERP_KERNEL::Exception(oss.str().c_str());
2246     }
2247 }
2248
2249 std::string MEDCoupling1DGTUMesh::simpleRepr() const
2250 {
2251   static const char msg0[]="No coordinates specified !";
2252   std::ostringstream ret;
2253   ret << "Single dynamic geometic type (" << _cm->getRepr() << ") unstructured mesh with name : \"" << getName() << "\"\n";
2254   ret << "Description of mesh : \"" << getDescription() << "\"\n";
2255   int tmpp1,tmpp2;
2256   double tt=getTime(tmpp1,tmpp2);
2257   ret << "Time attached to the mesh [unit] : " << tt << " [" << getTimeUnit() << "]\n";
2258   ret << "Iteration : " << tmpp1  << " Order : " << tmpp2 << "\n";
2259   ret << "Mesh dimension : " << getMeshDimension() << "\nSpace dimension : ";
2260   if(_coords!=0)
2261     {
2262       const int spaceDim=getSpaceDimension();
2263       ret << spaceDim << "\nInfo attached on space dimension : ";
2264       for(int i=0;i<spaceDim;i++)
2265         ret << "\"" << _coords->getInfoOnComponent(i) << "\" ";
2266       ret << "\n";
2267     }
2268   else
2269     ret << msg0 << "\n";
2270   ret << "Number of nodes : ";
2271   if(_coords!=0)
2272     ret << getNumberOfNodes() << "\n";
2273   else
2274     ret << msg0 << "\n";
2275   ret << "Number of cells : ";
2276   bool isOK=true;
2277   try { checkCoherency(); } catch(INTERP_KERNEL::Exception& e)
2278     {
2279       ret << "Nodal connectivity arrays are not set or badly set !\n";
2280       isOK=false;
2281     }
2282   if(isOK)
2283     ret << getNumberOfCells() << "\n";
2284   ret << "Cell type : " << _cm->getRepr() << "\n";
2285   return ret.str();
2286 }
2287
2288 std::string MEDCoupling1DGTUMesh::advancedRepr() const
2289 {
2290   std::ostringstream ret;
2291   ret << simpleRepr();
2292   ret << "\nCoordinates array : \n___________________\n\n";
2293   if(_coords)
2294     _coords->reprWithoutNameStream(ret);
2295   else
2296     ret << "No array set !\n";
2297   ret << "\n\nNodal Connectivity : \n____________________\n\n";
2298   //
2299   bool isOK=true;
2300   try { checkCoherency1(); } catch(INTERP_KERNEL::Exception& e)
2301     {
2302       ret << "Nodal connectivity arrays are not set or badly set !\n";
2303       isOK=false;
2304     }
2305   if(!isOK)
2306     return ret.str();
2307   int nbOfCells=getNumberOfCells();
2308   const int *ci=_conn_indx->begin(),*c=_conn->begin();
2309   for(int i=0;i<nbOfCells;i++,ci++)
2310     {
2311       ret << "Cell #" << i << " : ";
2312       std::copy(c+ci[0],c+ci[1],std::ostream_iterator<int>(ret," "));
2313       ret << "\n";
2314     }
2315   return ret.str();
2316 }
2317
2318 DataArrayDouble *MEDCoupling1DGTUMesh::computeIsoBarycenterOfNodesPerCell() const throw(INTERP_KERNEL::Exception)
2319 {
2320   MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> ret=DataArrayDouble::New();
2321   int spaceDim=getSpaceDimension();
2322   int nbOfCells=getNumberOfCells();//checkCoherency()
2323   int nbOfNodes=getNumberOfNodes();
2324   ret->alloc(nbOfCells,spaceDim);
2325   double *ptToFill=ret->getPointer();
2326   const double *coor=_coords->begin();
2327   const int *nodal=_conn->begin(),*nodali=_conn_indx->begin();
2328   nodal+=nodali[0];
2329   if(getCellModelEnum()!=INTERP_KERNEL::NORM_POLYHED)
2330     {
2331       for(int i=0;i<nbOfCells;i++,ptToFill+=spaceDim,nodali++)
2332         {
2333           std::fill(ptToFill,ptToFill+spaceDim,0.);
2334           if(nodali[0]<nodali[1])// >= to avoid division by 0.
2335             {
2336               for(int j=nodali[0];j<nodali[1];j++,nodal++)
2337                 {
2338                   if(*nodal>=0 && *nodal<nbOfNodes)
2339                     std::transform(coor+spaceDim*nodal[0],coor+spaceDim*(nodal[0]+1),ptToFill,ptToFill,std::plus<double>());
2340                   else
2341                     {
2342                       std::ostringstream oss; oss << "MEDCoupling1DGTUMesh::computeIsoBarycenterOfNodesPerCell : on cell #" << i << " presence of nodeId #" << *nodal << " should be in [0," <<   nbOfNodes << ") !";
2343                       throw INTERP_KERNEL::Exception(oss.str().c_str());
2344                     }
2345                   std::transform(ptToFill,ptToFill+spaceDim,ptToFill,std::bind2nd(std::multiplies<double>(),1./(nodali[1]-nodali[0])));
2346                 }
2347             }
2348           else
2349             {
2350               std::ostringstream oss; oss << "MEDCoupling1DGTUMesh::computeIsoBarycenterOfNodesPerCell : at cell #" << i << " the nodal index array is invalid !";
2351               throw INTERP_KERNEL::Exception(oss.str().c_str());
2352             }
2353         }
2354     }
2355   else
2356     {
2357       for(int i=0;i<nbOfCells;i++,ptToFill+=spaceDim,nodali++)
2358         {
2359           std::fill(ptToFill,ptToFill+spaceDim,0.);
2360           if(nodali[0]<nodali[1])// >= to avoid division by 0.
2361             {
2362               int nbOfNod=0;
2363               for(int j=nodali[0];j<nodali[1];j++,nodal++)
2364                 {
2365                   if(*nodal==-1) continue;
2366                   if(*nodal>=0 && *nodal<nbOfNodes)
2367                     {
2368                       std::transform(coor+spaceDim*nodal[0],coor+spaceDim*(nodal[0]+1),ptToFill,ptToFill,std::plus<double>());
2369                       nbOfNod++;
2370                     }
2371                   else
2372                     {
2373                       std::ostringstream oss; oss << "MEDCoupling1DGTUMesh::computeIsoBarycenterOfNodesPerCell (polyhedron) : on cell #" << i << " presence of nodeId #" << *nodal << " should be in [0," <<   nbOfNodes << ") !";
2374                       throw INTERP_KERNEL::Exception(oss.str().c_str());
2375                     }
2376                 }
2377               if(nbOfNod!=0)
2378                 std::transform(ptToFill,ptToFill+spaceDim,ptToFill,std::bind2nd(std::multiplies<double>(),1./nbOfNod));
2379               else
2380                 {
2381                   std::ostringstream oss; oss << "MEDCoupling1DGTUMesh::computeIsoBarycenterOfNodesPerCell (polyhedron) : no nodes in cell #" << i << " !";
2382                   throw INTERP_KERNEL::Exception(oss.str().c_str());
2383                 }
2384             }
2385           else
2386             {
2387               std::ostringstream oss; oss << "MEDCoupling1DGTUMesh::computeIsoBarycenterOfNodesPerCell (polyhedron)  : at cell #" << i << " the nodal index array is invalid !";
2388               throw INTERP_KERNEL::Exception(oss.str().c_str());
2389             }
2390         }
2391     }
2392   return ret.retn();
2393 }
2394
2395 void MEDCoupling1DGTUMesh::renumberCells(const int *old2NewBg, bool check) throw(INTERP_KERNEL::Exception)
2396 {
2397   int nbCells=getNumberOfCells();
2398   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> o2n=DataArrayInt::New();
2399   o2n->useArray(old2NewBg,false,C_DEALLOC,nbCells,1);
2400   if(check)
2401     o2n=o2n->checkAndPreparePermutation();
2402   //
2403   const int *o2nPtr=o2n->getPointer();
2404   const int *conn=_conn->begin(),*conni=_conn_indx->begin();
2405   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> newConn=DataArrayInt::New();
2406   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> newConnI=DataArrayInt::New();
2407   newConn->alloc(_conn->getNumberOfTuples(),1); newConnI->alloc(nbCells,1);
2408   newConn->copyStringInfoFrom(*_conn); newConnI->copyStringInfoFrom(*_conn_indx);
2409   //
2410   int *newC=newConn->getPointer(),*newCI=newConnI->getPointer();
2411   for(int i=0;i<nbCells;i++)
2412     {
2413       int newPos=o2nPtr[i];
2414       int sz=conni[i+1]-conni[i];
2415       if(sz>=0)
2416         newCI[newPos]=sz;
2417       else
2418         {
2419           std::ostringstream oss; oss << "MEDCoupling1DGTUMesh::renumberCells : the index nodal array is invalid for cell #" << i << " !";
2420           throw INTERP_KERNEL::Exception(oss.str().c_str());
2421         }
2422     }
2423   newConnI->computeOffsets2(); newCI=newConnI->getPointer();
2424   //
2425   for(int i=0;i<nbCells;i++,conni++)
2426     {
2427       int sz=conni[1]-conni[0];
2428       int newp=o2nPtr[i];
2429       std::copy(conn+conni[0],conn+conni[1],newC+newCI[newp]);
2430     }
2431   _conn=newConn;
2432   _conn_indx=newConnI;
2433 }
2434
2435 MEDCouplingMesh *MEDCoupling1DGTUMesh::mergeMyselfWith(const MEDCouplingMesh *other) const
2436 {
2437   if(other->getType()!=SINGLE_DYNAMIC_GEO_TYPE_UNSTRUCTURED)
2438     throw INTERP_KERNEL::Exception("Merge of umesh only available with umesh single dynamic geo type each other !");
2439   const MEDCoupling1DGTUMesh *otherC=static_cast<const MEDCoupling1DGTUMesh *>(other);
2440   return Merge1DGTUMeshes(this,otherC);
2441 }
2442
2443 MEDCouplingUMesh *MEDCoupling1DGTUMesh::buildUnstructured() const throw(INTERP_KERNEL::Exception)
2444 {
2445   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> ret=MEDCouplingUMesh::New(getName().c_str(),getMeshDimension());
2446   ret->setCoords(getCoords());
2447   const int *nodalConn=_conn->begin(),*nodalConnI=_conn_indx->begin();
2448   int nbCells=getNumberOfCells();//checkCoherency
2449   int geoType=(int)getCellModelEnum();
2450   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> c=DataArrayInt::New(); c->alloc(nbCells+_conn->getNumberOfTuples(),1);
2451   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> cI=DataArrayInt::New(); cI->alloc(nbCells+1);
2452   int *cPtr=c->getPointer(),*ciPtr=cI->getPointer();
2453   ciPtr[0]=0;
2454   for(int i=0;i<nbCells;i++,ciPtr++)
2455     {
2456       int sz=nodalConnI[i+1]-nodalConnI[i];
2457       if(sz>=0)
2458         {
2459           *cPtr++=geoType;
2460           cPtr=std::copy(nodalConn+nodalConnI[i],nodalConn+nodalConnI[i+1],cPtr);
2461           ciPtr[1]=ciPtr[0]+sz+1;
2462         }
2463       else
2464         {
2465           std::ostringstream oss; oss << "MEDCoupling1DGTUMesh::buildUnstructured : Invalid for nodal index for cell #" << i << " !";
2466           throw INTERP_KERNEL::Exception(oss.str().c_str());
2467         }
2468     }
2469   ret->setConnectivity(c,cI,true);
2470   return ret.retn();
2471 }
2472
2473 /*!
2474  * Do nothing for the moment, because there is no policy that allows to split polygons, polyhedrons ... into simplexes
2475  */
2476 DataArrayInt *MEDCoupling1DGTUMesh::simplexize(int policy) throw(INTERP_KERNEL::Exception)
2477 {
2478   int nbOfCells=getNumberOfCells();
2479   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret=DataArrayInt::New();
2480   ret->alloc(nbOfCells,1);
2481   ret->iota(0);
2482   return ret.retn();
2483 }
2484
2485 void MEDCoupling1DGTUMesh::reprQuickOverview(std::ostream& stream) const throw(INTERP_KERNEL::Exception)
2486 {
2487   stream << "MEDCoupling1DGTUMesh C++ instance at " << this << ". Type=" << _cm->getRepr() << ". Name : \"" << getName() << "\".";
2488   stream << " Mesh dimension : " << getMeshDimension() << ".";
2489   if(!_coords)
2490     { stream << " No coordinates set !"; return ; }
2491   if(!_coords->isAllocated())
2492     { stream << " Coordinates set but not allocated !"; return ; }
2493   stream << " Space dimension : " << _coords->getNumberOfComponents() << "." << std::endl;
2494   stream << "Number of nodes : " << _coords->getNumberOfTuples() << ".";
2495   bool isOK=true;
2496   try { checkCoherency(); } catch(INTERP_KERNEL::Exception& e)
2497     {
2498       stream << std::endl << "Nodal connectivity NOT set properly !\n";
2499       isOK=false;
2500     }
2501   if(isOK)
2502     stream << std::endl << "Number of cells : " << getNumberOfCells() << ".";
2503 }
2504
2505 void MEDCoupling1DGTUMesh::shallowCopyConnectivityFrom(const MEDCouplingPointSet *other) throw(INTERP_KERNEL::Exception)
2506 {
2507   if(!other)
2508     throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::shallowCopyConnectivityFrom : input pointer is null !");
2509   const MEDCoupling1DGTUMesh *otherC=dynamic_cast<const MEDCoupling1DGTUMesh *>(other);
2510   if(!otherC)
2511     throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::shallowCopyConnectivityFrom : input pointer is not an MEDCoupling1DGTUMesh instance !");
2512   setNodalConnectivity(otherC->getNodalConnectivity(),otherC->getNodalConnectivityIndex());
2513 }
2514
2515 MEDCouplingPointSet *MEDCoupling1DGTUMesh::mergeMyselfWithOnSameCoords(const MEDCouplingPointSet *other) const
2516 {
2517   if(!other)
2518     throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::mergeMyselfWithOnSameCoords : input other is null !");
2519   const MEDCoupling1DGTUMesh *otherC=dynamic_cast<const MEDCoupling1DGTUMesh *>(other);
2520   if(!otherC)
2521     throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::mergeMyselfWithOnSameCoords : the input other mesh is not of type single statuc geo type unstructured !");
2522   std::vector<const MEDCoupling1DGTUMesh *> ms(2);
2523   ms[0]=this;
2524   ms[1]=otherC;
2525   return Merge1DGTUMeshesOnSameCoords(ms);
2526 }
2527
2528 MEDCouplingPointSet *MEDCoupling1DGTUMesh::buildPartOfMySelfKeepCoords(const int *begin, const int *end) const
2529 {
2530   checkCoherency();
2531   MEDCouplingAutoRefCountObjectPtr<MEDCoupling1DGTUMesh> ret(new MEDCoupling1DGTUMesh(getName().c_str(),*_cm));
2532   ret->setCoords(_coords);
2533   DataArrayInt *c=0,*ci=0;
2534   MEDCouplingUMesh::ExtractFromIndexedArrays(begin,end,_conn,_conn_indx,c,ci);
2535   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> cSafe(c),ciSafe(ci);
2536   ret->setNodalConnectivity(c,ci);
2537   return ret.retn();
2538 }
2539
2540 MEDCouplingPointSet *MEDCoupling1DGTUMesh::buildPartOfMySelfKeepCoords2(int start, int end, int step) const
2541 {
2542   checkCoherency();
2543   MEDCouplingAutoRefCountObjectPtr<MEDCoupling1DGTUMesh> ret(new MEDCoupling1DGTUMesh(getName().c_str(),*_cm));
2544   ret->setCoords(_coords);
2545   DataArrayInt *c=0,*ci=0;
2546   MEDCouplingUMesh::ExtractFromIndexedArrays2(start,end,step,_conn,_conn_indx,c,ci);
2547   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> cSafe(c),ciSafe(ci);
2548   ret->setNodalConnectivity(c,ci);
2549   return ret.retn();
2550 }
2551
2552 void MEDCoupling1DGTUMesh::computeNodeIdsAlg(std::vector<bool>& nodeIdsInUse) const throw(INTERP_KERNEL::Exception)
2553 {
2554   int sz((int)nodeIdsInUse.size());
2555   int nbCells(getNumberOfCells());
2556   const int *w(_conn->begin()),*wi(_conn_indx->begin());
2557   for(int i=0;i<nbCells;i++,wi++)
2558     for(const int *pt=w+wi[0];pt!=w+wi[1];pt++)
2559       if(*pt!=-1)
2560         {
2561           if(*pt>=0 && *pt<sz)
2562             nodeIdsInUse[*pt]=true;
2563           else
2564             {
2565               std::ostringstream oss; oss << "MEDCoupling1DGTUMesh::computeNodeIdsAlg : At cell #" << i << " presence of node id #" << *pt << " should be in [0," << sz << ") !";
2566               throw INTERP_KERNEL::Exception(oss.str().c_str());
2567             }
2568         }
2569 }
2570
2571 void MEDCoupling1DGTUMesh::getReverseNodalConnectivity(DataArrayInt *revNodal, DataArrayInt *revNodalIndx) const throw(INTERP_KERNEL::Exception)
2572 {
2573   checkFullyDefined();
2574   int nbOfNodes=getNumberOfNodes();
2575   int *revNodalIndxPtr=(int *)malloc((nbOfNodes+1)*sizeof(int));
2576   revNodalIndx->useArray(revNodalIndxPtr,true,C_DEALLOC,nbOfNodes+1,1);
2577   std::fill(revNodalIndxPtr,revNodalIndxPtr+nbOfNodes+1,0);
2578   const int *conn=_conn->begin(),*conni=_conn_indx->begin();
2579   int nbOfCells=getNumberOfCells();
2580   int nbOfEltsInRevNodal=0;
2581   for(int eltId=0;eltId<nbOfCells;eltId++)
2582     {
2583       int nbOfNodesPerCell=conni[eltId+1]-conni[eltId];
2584       if(nbOfNodesPerCell>=0)
2585         {
2586           for(int j=0;j<nbOfNodesPerCell;j++)
2587             {
2588               int nodeId=conn[conni[eltId]+j];
2589               if(nodeId==-1) continue;
2590               if(nodeId>=0 && nodeId<nbOfNodes)
2591                 {
2592                   nbOfEltsInRevNodal++;
2593                   revNodalIndxPtr[nodeId+1]++;
2594                 }
2595               else
2596                 {
2597                   std::ostringstream oss; oss << "MEDCoupling1DGTUMesh::getReverseNodalConnectivity : At cell #" << eltId << " presence of nodeId #" << conn[0] << " should be in [0," << nbOfNodes << ") !";
2598                   throw INTERP_KERNEL::Exception(oss.str().c_str());
2599                 }
2600             }
2601         }
2602       else
2603         {
2604           std::ostringstream oss; oss << "MEDCoupling1DGTUMesh::getReverseNodalConnectivity : At cell #" << eltId << "nodal connectivity is invalid !";
2605           throw INTERP_KERNEL::Exception(oss.str().c_str());
2606         }
2607     }
2608   std::transform(revNodalIndxPtr+1,revNodalIndxPtr+nbOfNodes+1,revNodalIndxPtr,revNodalIndxPtr+1,std::plus<int>());
2609   conn=_conn->begin();
2610   int *revNodalPtr=(int *)malloc((nbOfEltsInRevNodal)*sizeof(int));
2611   revNodal->useArray(revNodalPtr,true,C_DEALLOC,nbOfEltsInRevNodal,1);
2612   std::fill(revNodalPtr,revNodalPtr+nbOfEltsInRevNodal,-1);
2613   for(int eltId=0;eltId<nbOfCells;eltId++)
2614     {
2615       int nbOfNodesPerCell=conni[eltId+1]-conni[eltId];
2616       for(int j=0;j<nbOfNodesPerCell;j++)
2617         {
2618           int nodeId=conn[conni[eltId]+j];
2619           if(nodeId!=-1)
2620             *std::find_if(revNodalPtr+revNodalIndxPtr[nodeId],revNodalPtr+revNodalIndxPtr[nodeId+1],std::bind2nd(std::equal_to<int>(),-1))=eltId;
2621         }
2622     }
2623 }
2624
2625 void MEDCoupling1DGTUMesh::checkFullyDefined() const throw(INTERP_KERNEL::Exception)
2626 {
2627   if(!((const DataArrayInt *)_conn) || !((const DataArrayInt *)_conn_indx) || !((const DataArrayDouble *)_coords))
2628     throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::checkFullyDefined : part of this is not fully defined.");
2629 }
2630
2631 bool MEDCoupling1DGTUMesh::isEmptyMesh(const std::vector<int>& tinyInfo) const
2632 {
2633   throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::isEmptyMesh : not implemented yet !");
2634 }
2635
2636 void MEDCoupling1DGTUMesh::getTinySerializationInformation(std::vector<double>& tinyInfoD, std::vector<int>& tinyInfo, std::vector<std::string>& littleStrings) const
2637 {
2638   int it,order;
2639   double time=getTime(it,order);
2640   tinyInfo.clear(); tinyInfoD.clear(); littleStrings.clear();
2641   //
2642   littleStrings.push_back(getName());
2643   littleStrings.push_back(getDescription());
2644   littleStrings.push_back(getTimeUnit());
2645   //
2646   std::vector<std::string> littleStrings2,littleStrings3,littleStrings4;
2647   if((const DataArrayDouble *)_coords)
2648     _coords->getTinySerializationStrInformation(littleStrings2);
2649   if((const DataArrayInt *)_conn)
2650     _conn->getTinySerializationStrInformation(littleStrings3);
2651   if((const DataArrayInt *)_conn_indx)
2652     _conn_indx->getTinySerializationStrInformation(littleStrings4);
2653   int sz0((int)littleStrings2.size()),sz1((int)littleStrings3.size()),sz2((int)littleStrings4.size());
2654   littleStrings.insert(littleStrings.end(),littleStrings2.begin(),littleStrings2.end());
2655   littleStrings.insert(littleStrings.end(),littleStrings3.begin(),littleStrings3.end());
2656   littleStrings.insert(littleStrings.end(),littleStrings4.begin(),littleStrings4.end());
2657   //
2658   tinyInfo.push_back(getCellModelEnum());
2659   tinyInfo.push_back(it);
2660   tinyInfo.push_back(order);
2661   std::vector<int> tinyInfo2,tinyInfo3,tinyInfo4;
2662   if((const DataArrayDouble *)_coords)
2663     _coords->getTinySerializationIntInformation(tinyInfo2);
2664   if((const DataArrayInt *)_conn)
2665     _conn->getTinySerializationIntInformation(tinyInfo3);
2666   if((const DataArrayInt *)_conn_indx)
2667     _conn_indx->getTinySerializationIntInformation(tinyInfo4);
2668   int sz3((int)tinyInfo2.size()),sz4((int)tinyInfo3.size()),sz5((int)tinyInfo4.size());
2669   tinyInfo.push_back(sz0); tinyInfo.push_back(sz1); tinyInfo.push_back(sz2); tinyInfo.push_back(sz3); tinyInfo.push_back(sz4);  tinyInfo.push_back(sz5);
2670   tinyInfo.insert(tinyInfo.end(),tinyInfo2.begin(),tinyInfo2.end());
2671   tinyInfo.insert(tinyInfo.end(),tinyInfo3.begin(),tinyInfo3.end());
2672   tinyInfo.insert(tinyInfo.end(),tinyInfo4.begin(),tinyInfo4.end());
2673   //
2674   tinyInfoD.push_back(time);
2675 }
2676
2677 void MEDCoupling1DGTUMesh::resizeForUnserialization(const std::vector<int>& tinyInfo, DataArrayInt *a1, DataArrayDouble *a2, std::vector<std::string>& littleStrings) const
2678 {
2679   std::vector<int> tinyInfo2(tinyInfo.begin()+9,tinyInfo.begin()+9+tinyInfo[6]);
2680   std::vector<int> tinyInfo1(tinyInfo.begin()+9+tinyInfo[6],tinyInfo.begin()+9+tinyInfo[6]+tinyInfo[7]);
2681   std::vector<int> tinyInfo12(tinyInfo.begin()+9+tinyInfo[6]+tinyInfo[7],tinyInfo.begin()+9+tinyInfo[6]+tinyInfo[7]+tinyInfo[8]);
2682   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> p1(DataArrayInt::New()); p1->resizeForUnserialization(tinyInfo1);
2683   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> p2(DataArrayInt::New()); p2->resizeForUnserialization(tinyInfo12);
2684   std::vector<const DataArrayInt *> v(2); v[0]=p1; v[1]=p2;
2685   p2=DataArrayInt::Aggregate(v);
2686   a2->resizeForUnserialization(tinyInfo2);
2687   a1->alloc(p2->getNbOfElems(),1);
2688 }
2689
2690 void MEDCoupling1DGTUMesh::serialize(DataArrayInt *&a1, DataArrayDouble *&a2) const
2691 {
2692   int sz(0);
2693   if((const DataArrayInt *)_conn)
2694     if(_conn->isAllocated())
2695       sz=_conn->getNbOfElems();
2696   if((const DataArrayInt *)_conn_indx)
2697     if(_conn_indx->isAllocated())
2698       sz+=_conn_indx->getNbOfElems();
2699   a1=DataArrayInt::New();
2700   a1->alloc(sz,1);
2701   int *work(a1->getPointer());
2702   if(sz!=0 && (const DataArrayInt *)_conn)
2703     work=std::copy(_conn->begin(),_conn->end(),a1->getPointer());
2704   if(sz!=0 && (const DataArrayInt *)_conn_indx)
2705     std::copy(_conn_indx->begin(),_conn_indx->end(),work);
2706   sz=0;
2707   if((const DataArrayDouble *)_coords)
2708     if(_coords->isAllocated())
2709       sz=_coords->getNbOfElems();
2710   a2=DataArrayDouble::New();
2711   a2->alloc(sz,1);
2712   if(sz!=0 && (const DataArrayDouble *)_coords)
2713     std::copy(_coords->begin(),_coords->end(),a2->getPointer());
2714 }
2715
2716 void MEDCoupling1DGTUMesh::unserialization(const std::vector<double>& tinyInfoD, const std::vector<int>& tinyInfo, const DataArrayInt *a1, DataArrayDouble *a2,
2717                                            const std::vector<std::string>& littleStrings)
2718 {
2719   INTERP_KERNEL::NormalizedCellType gt((INTERP_KERNEL::NormalizedCellType)tinyInfo[0]);
2720   _cm=&INTERP_KERNEL::CellModel::GetCellModel(gt);
2721   setName(littleStrings[0].c_str());
2722   setDescription(littleStrings[1].c_str());
2723   setTimeUnit(littleStrings[2].c_str());
2724   setTime(tinyInfoD[0],tinyInfo[1],tinyInfo[2]);
2725   int sz0(tinyInfo[3]),sz1(tinyInfo[4]),sz2(tinyInfo[5]),sz3(tinyInfo[6]),sz4(tinyInfo[7]),sz5(tinyInfo[8]);
2726   //
2727   _coords=DataArrayDouble::New();
2728   std::vector<int> tinyInfo2(tinyInfo.begin()+9,tinyInfo.begin()+9+sz3);
2729   _coords->resizeForUnserialization(tinyInfo2);
2730   std::copy(a2->begin(),a2->end(),_coords->getPointer());
2731   _conn=DataArrayInt::New();
2732   std::vector<int> tinyInfo3(tinyInfo.begin()+9+sz3,tinyInfo.begin()+9+sz3+sz4);
2733   _conn->resizeForUnserialization(tinyInfo3);
2734   std::copy(a1->begin(),a1->begin()+_conn->getNbOfElems(),_conn->getPointer());
2735   _conn_indx=DataArrayInt::New();
2736   std::vector<int> tinyInfo4(tinyInfo.begin()+9+sz3+sz4,tinyInfo.begin()+9+sz3+sz4+sz5);
2737   _conn_indx->resizeForUnserialization(tinyInfo4);
2738   std::copy(a1->begin()+_conn->getNbOfElems(),a1->end(),_conn_indx->getPointer());
2739   std::vector<std::string> littleStrings2(littleStrings.begin()+3,littleStrings.begin()+3+sz0);
2740   _coords->finishUnserialization(tinyInfo2,littleStrings2);
2741   std::vector<std::string> littleStrings3(littleStrings.begin()+3+sz0,littleStrings.begin()+3+sz0+sz1);
2742   _conn->finishUnserialization(tinyInfo3,littleStrings3);
2743   std::vector<std::string> littleStrings4(littleStrings.begin()+3+sz0+sz1,littleStrings.begin()+3+sz0+sz1+sz2);
2744   _conn_indx->finishUnserialization(tinyInfo4,littleStrings4);
2745 }
2746
2747 /*!
2748  * Finds nodes not used in any cell and returns an array giving a new id to every node
2749  * by excluding the unused nodes, for which the array holds -1. The result array is
2750  * a mapping in "Old to New" mode.
2751  *  \param [out] nbrOfNodesInUse - number of node ids present in the nodal connectivity.
2752  *  \return DataArrayInt * - a new instance of DataArrayInt. Its length is \a
2753  *          this->getNumberOfNodes(). It holds for each node of \a this mesh either -1
2754  *          if the node is unused or a new id else. The caller is to delete this
2755  *          array using decrRef() as it is no more needed.
2756  *  \throw If the coordinates array is not set.
2757  *  \throw If the nodal connectivity of cells is not defined.
2758  *  \throw If the nodal connectivity includes an invalid id.
2759  */
2760 DataArrayInt *MEDCoupling1DGTUMesh::getNodeIdsInUse(int& nbrOfNodesInUse) const throw(INTERP_KERNEL::Exception)
2761 {
2762   nbrOfNodesInUse=-1;
2763   int nbOfNodes=getNumberOfNodes();
2764   int nbOfCells=getNumberOfCells();//checkCoherency
2765   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret=DataArrayInt::New();
2766   ret->alloc(nbOfNodes,1);
2767   int *traducer=ret->getPointer();
2768   std::fill(traducer,traducer+nbOfNodes,-1);
2769   const int *conn=_conn->begin(),*conni(_conn_indx->begin());
2770   for(int i=0;i<nbOfCells;i++,conni++)
2771     {
2772       int nbNodesPerCell=conni[1]-conni[0];
2773       for(int j=0;j<nbNodesPerCell;j++)
2774         {
2775           int nodeId=conn[conni[0]+j];
2776           if(nodeId==-1) continue;
2777           if(nodeId>=0 && nodeId<nbOfNodes)
2778             traducer[nodeId]=1;
2779           else
2780             {
2781               std::ostringstream oss; oss << "MEDCoupling1DGTUMesh::getNodeIdsInUse : In cell #" << i  << " presence of node id " <<  nodeId << " not in [0," << nbOfNodes << ") !";
2782               throw INTERP_KERNEL::Exception(oss.str().c_str());
2783             }
2784         }
2785     }
2786   nbrOfNodesInUse=(int)std::count(traducer,traducer+nbOfNodes,1);
2787   std::transform(traducer,traducer+nbOfNodes,traducer,MEDCouplingAccVisit());
2788   return ret.retn();
2789 }
2790
2791 /*!
2792  * Changes ids of nodes within the nodal connectivity arrays according to a permutation
2793  * array in "Old to New" mode. The node coordinates array is \b not changed by this method.
2794  * This method is a generalization of shiftNodeNumbersInConn().
2795  *  \warning This method performs no check of validity of new ids. **Use it with care !**
2796  *  \param [in] newNodeNumbersO2N - a permutation array, of length \a
2797  *         this->getNumberOfNodes(), in "Old to New" mode.
2798  *         See \ref MEDCouplingArrayRenumbering for more info on renumbering modes.
2799  *  \throw If the nodal connectivity of cells is not defined.
2800  */
2801 void MEDCoupling1DGTUMesh::renumberNodesInConn(const int *newNodeNumbersO2N)
2802 {
2803   getNumberOfCells();//only to check that all is well defined.
2804   //
2805   int nbElemsIn=getNumberOfNodes();
2806   int nbOfTuples=_conn->getNumberOfTuples();
2807   int *pt=_conn->getPointer();
2808   for(int i=0;i<nbOfTuples;i++,pt++)
2809     {
2810       if(*pt==-1) continue;
2811       if(*pt>=0 && *pt<nbElemsIn)
2812         *pt=newNodeNumbersO2N[*pt];
2813       else
2814         {
2815           std::ostringstream oss; oss << "MEDCoupling1DGTUMesh::renumberNodesInConn : error on tuple #" << i << " value is " << *pt << " and indirectionnal array as a size equal to " << nbElemsIn;
2816           throw INTERP_KERNEL::Exception(oss.str().c_str());
2817         }
2818     }
2819   _conn->declareAsNew();
2820   //
2821   updateTime();
2822 }
2823
2824 /*!
2825  * Keeps from \a this only cells which constituing point id are in the ids specified by [\a begin,\a end).
2826  * The resulting cell ids are stored at the end of the 'cellIdsKept' parameter.
2827  * Parameter \a fullyIn specifies if a cell that has part of its nodes in ids array is kept or not.
2828  * If \a fullyIn is true only cells whose ids are \b fully contained in [\a begin,\a end) tab will be kept.
2829  *
2830  * \param [in] begin input start of array of node ids.
2831  * \param [in] end input end of array of node ids.
2832  * \param [in] fullyIn input that specifies if all node ids must be in [\a begin,\a end) array to consider cell to be in.
2833  * \param [in,out] cellIdsKeptArr array where all candidate cell ids are put at the end.
2834  */
2835 void MEDCoupling1DGTUMesh::fillCellIdsToKeepFromNodeIds(const int *begin, const int *end, bool fullyIn, DataArrayInt *&cellIdsKeptArr) const
2836 {
2837   int nbOfCells=getNumberOfCells();
2838   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> cellIdsKept=DataArrayInt::New(); cellIdsKept->alloc(0,1);
2839   int tmp=-1;
2840   int sz=_conn->getMaxValue(tmp); sz=std::max(sz,0)+1;
2841   std::vector<bool> fastFinder(sz,false);
2842   for(const int *work=begin;work!=end;work++)
2843     if(*work>=0 && *work<sz)
2844       fastFinder[*work]=true;
2845   const int *conn=_conn->begin(),*conni=_conn_indx->begin();
2846   for(int i=0;i<nbOfCells;i++,conni++)
2847     {
2848       int ref=0,nbOfHit=0;
2849       int nbNodesPerCell=conni[1]-conni[0];
2850       if(nbNodesPerCell>=0)
2851         {
2852           for(int j=0;j<nbNodesPerCell;j++)
2853             {
2854               int nodeId=conn[conni[0]+j];
2855               if(nodeId>=0)
2856                 {
2857                   ref++;
2858                   if(fastFinder[nodeId])
2859                     nbOfHit++;
2860                 }
2861             }
2862         }
2863       else
2864         {
2865           std::ostringstream oss; oss << "MEDCoupling1DGTUMesh::fillCellIdsToKeepFromNodeIds : invalid index array for cell #" << i << " !";
2866           throw INTERP_KERNEL::Exception(oss.str().c_str());
2867         }
2868       if((ref==nbOfHit && fullyIn) || (nbOfHit!=0 && !fullyIn))
2869         cellIdsKept->pushBackSilent(i);
2870     }
2871   cellIdsKeptArr=cellIdsKept.retn();
2872 }
2873
2874 void MEDCoupling1DGTUMesh::allocateCells(int nbOfCells) throw(INTERP_KERNEL::Exception)
2875 {
2876   if(nbOfCells<0)
2877     throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::allocateCells : the input number of cells should be >= 0 !");
2878   _conn=DataArrayInt::New();
2879   _conn->reserve(nbOfCells*3);
2880   _conn_indx=DataArrayInt::New();
2881   _conn_indx->reserve(nbOfCells+1); _conn_indx->pushBackSilent(0);
2882   declareAsNew();
2883 }
2884
2885 /*!
2886  * Appends at the end of \a this a cell having nodal connectivity array defined in [ \a nodalConnOfCellBg, \a nodalConnOfCellEnd ).
2887  *
2888  * \param [in] nodalConnOfCellBg - the begin (included) of nodal connectivity of the cell to add.
2889  * \param [in] nodalConnOfCellEnd - the end (excluded) of nodal connectivity of the cell to add.
2890  * \throw If the length of the input nodal connectivity array of the cell to add is not equal to number of nodes per cell relative to the unique geometric type
2891  *        attached to \a this.
2892  * \thow If the nodal connectivity array in \a this is null (call MEDCoupling1SGTUMesh::allocateCells before).
2893  */
2894 void MEDCoupling1DGTUMesh::insertNextCell(const int *nodalConnOfCellBg, const int *nodalConnOfCellEnd) throw(INTERP_KERNEL::Exception)
2895 {
2896   int sz=(int)std::distance(nodalConnOfCellBg,nodalConnOfCellEnd);
2897   DataArrayInt *c(_conn),*c2(_conn_indx);
2898   if(c && c2)
2899     {
2900       int pos=c2->back();
2901       if(pos==c->getNumberOfTuples())
2902         {
2903           c->pushBackValsSilent(nodalConnOfCellBg,nodalConnOfCellEnd);
2904           c2->pushBackSilent(pos+sz);
2905         }
2906       else
2907         {
2908           std::ostringstream oss; oss << "MEDCoupling1DGTUMesh::insertNextCell : The nodal index array (end=" << pos << ") mismatches with nodal array (length=" << c->getNumberOfTuples() << ") !";
2909           throw INTERP_KERNEL::Exception(oss.str().c_str());
2910         }
2911     }
2912   else
2913     throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::insertNextCell : nodal connectivity array is null ! Call MEDCoupling1DGTUMesh::allocateCells before !");
2914 }
2915
2916 void MEDCoupling1DGTUMesh::setNodalConnectivity(DataArrayInt *nodalConn, DataArrayInt *nodalConnIndex) throw(INTERP_KERNEL::Exception)
2917 {
2918   if(nodalConn)
2919     nodalConn->incrRef();
2920   _conn=nodalConn;
2921   if(nodalConnIndex)
2922     nodalConnIndex->incrRef();
2923   _conn_indx=nodalConnIndex;
2924   declareAsNew();
2925 }
2926
2927 /*!
2928  * \return DataArrayInt * - the internal reference to the nodal connectivity. The caller is not reponsible to deallocate it.
2929  */
2930 DataArrayInt *MEDCoupling1DGTUMesh::getNodalConnectivity() const throw(INTERP_KERNEL::Exception)
2931 {
2932   const DataArrayInt *ret(_conn);
2933   return const_cast<DataArrayInt *>(ret);
2934 }
2935
2936 /*!
2937  * \return DataArrayInt * - the internal reference to the nodal connectivity index. The caller is not reponsible to deallocate it.
2938  */
2939 DataArrayInt *MEDCoupling1DGTUMesh::getNodalConnectivityIndex() const throw(INTERP_KERNEL::Exception)
2940 {
2941   const DataArrayInt *ret(_conn_indx);
2942   return const_cast<DataArrayInt *>(ret);
2943 }
2944
2945 /*!
2946  * See the definition of the nodal connectivity pack \ref MEDCoupling1DGTUMesh::isPacked "here".
2947  * This method tries to build a new instance geometrically equivalent to \a this, by limiting at most the number of new object (nodal connectivity).
2948  * Geometrically the returned mesh is equal to \a this. So if \a this is already packed, the return value is a shallow copy of \a this.
2949  *
2950  * Whatever the status of pack of \a this, the coordinates array of the returned newly created instance is the same than those in \a this.
2951  *
2952  * \param [out] isShallowCpyOfNodalConnn - tells if the returned instance share the same pair of nodal connectivity arrays (true) or if nodal
2953  *              connectivity arrays are different (false)
2954  * \return a new object to be managed by the caller.
2955  *
2956  * \sa MEDCoupling1DGTUMesh::retrievePackedNodalConnectivity, MEDCoupling1DGTUMesh::isPacked
2957  */
2958 MEDCoupling1DGTUMesh *MEDCoupling1DGTUMesh::copyWithNodalConnectivityPacked(bool& isShallowCpyOfNodalConnn) const throw(INTERP_KERNEL::Exception)
2959 {
2960   MEDCouplingAutoRefCountObjectPtr<MEDCoupling1DGTUMesh> ret(new MEDCoupling1DGTUMesh(getName().c_str(),*_cm));
2961   DataArrayInt *nc=0,*nci=0;
2962   isShallowCpyOfNodalConnn=retrievePackedNodalConnectivity(nc,nci);
2963   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ncs(nc),ncis(nci);
2964   ret->_conn=ncs; ret->_conn_indx=ncis;
2965   ret->setCoords(getCoords());
2966   return ret.retn();
2967 }
2968
2969 /*!
2970  * This method allows to compute, if needed, the packed nodal connectivity pair.
2971  * Indeed, it is possible to store in \a this a nodal connectivity array bigger than ranges convered by nodal connectivity index array.
2972  * It is typically the case when nodalConnIndx starts with an id greater than 0, and finishes with id less than number of tuples in \c this->_conn.
2973  *
2974  * If \a this looks packed (the front of nodal connectivity index equal to 0 and back of connectivity index equal to number of tuple of nodal connectivity array)
2975  * true will be returned and respectively \a this->_conn and \a this->_conn_indx (with ref counter incremented). This is the classical case.
2976  *
2977  * If nodal connectivity index points to a subpart of nodal connectivity index the packed pair of arrays will be computed (new objects) and returned and false
2978  * will be returned.
2979  *
2980  * This method return 3 elements.
2981  * \param [out] nodalConn - a pointer that can be equal to \a this->_conn if true is returned (general case). Whatever the value of return parameter
2982  *                          this pointer can be seen as a new object, that is to managed by the caller.
2983  * \param [out] nodalConnIndx - a pointer that can be equal to \a this->_conn_indx if true is returned (general case). Whatever the value of return parameter
2984  *                              this pointer can be seen as a new object, that is to managed by the caller.
2985  * \return bool - an indication of the content of the 2 output parameters. If true, \a this looks packed (general case), if true, \a this is not packed then
2986  * output parameters are newly created objects.
2987  *
2988  * \throw if \a this does not pass MEDCoupling1DGTUMesh::checkCoherency test
2989  */
2990 bool MEDCoupling1DGTUMesh::retrievePackedNodalConnectivity(DataArrayInt *&nodalConn, DataArrayInt *&nodalConnIndx) const throw(INTERP_KERNEL::Exception)
2991 {
2992   if(isPacked())//performs the checkCoherency
2993     {
2994       const DataArrayInt *c0(_conn),*c1(_conn_indx);
2995       nodalConn=const_cast<DataArrayInt *>(c0); nodalConnIndx=const_cast<DataArrayInt *>(c1);
2996       nodalConn->incrRef(); nodalConnIndx->incrRef();
2997       return true;
2998     }
2999   int bg=_conn_indx->front(),end=_conn_indx->back();
3000   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> nc(_conn->selectByTupleId2(bg,end,1));
3001   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> nci(_conn_indx->deepCpy());
3002   nci->applyLin(1,-bg);
3003   nodalConn=nc.retn(); nodalConnIndx=nci.retn();
3004   return false;
3005 }
3006
3007 /*
3008  * If \a this looks packed (the front of nodal connectivity index equal to 0 and back of connectivity index equal to number of tuple of nodal connectivity array)
3009  * true will be returned and respectively \a this->_conn and \a this->_conn_indx (with ref counter incremented). This is the classical case.
3010  * If nodal connectivity index points to a subpart of nodal connectivity index false will be returned.
3011  * \return bool - true if \a this looks packed, false is not.
3012  *
3013  * \throw if \a this does not pass MEDCoupling1DGTUMesh::checkCoherency test
3014  */
3015 bool MEDCoupling1DGTUMesh::isPacked() const throw(INTERP_KERNEL::Exception)
3016 {
3017   checkCoherency();
3018   return _conn_indx->front()==0 && _conn_indx->back()==_conn->getNumberOfTuples();
3019 }
3020
3021 MEDCoupling1DGTUMesh *MEDCoupling1DGTUMesh::Merge1DGTUMeshes(const MEDCoupling1DGTUMesh *mesh1, const MEDCoupling1DGTUMesh *mesh2) throw(INTERP_KERNEL::Exception)
3022 {
3023   std::vector<const MEDCoupling1DGTUMesh *> tmp(2);
3024   tmp[0]=const_cast<MEDCoupling1DGTUMesh *>(mesh1); tmp[1]=const_cast<MEDCoupling1DGTUMesh *>(mesh2);
3025   return Merge1DGTUMeshes(tmp);
3026 }
3027
3028 MEDCoupling1DGTUMesh *MEDCoupling1DGTUMesh::Merge1DGTUMeshes(std::vector<const MEDCoupling1DGTUMesh *>& a) throw(INTERP_KERNEL::Exception)
3029 {
3030   std::size_t sz=a.size();
3031   if(sz==0)
3032     return Merge1DGTUMeshesLL(a);
3033   for(std::size_t ii=0;ii<sz;ii++)
3034     if(!a[ii])
3035       {
3036         std::ostringstream oss; oss << "MEDCoupling1DGTUMesh::Merge1DGTUMeshes : item #" << ii << " in input array of size "<< sz << " is empty !";
3037         throw INTERP_KERNEL::Exception(oss.str().c_str());
3038       }
3039   const INTERP_KERNEL::CellModel *cm=&(a[0]->getCellModel());
3040   for(std::size_t ii=0;ii<sz;ii++)
3041     if(&(a[ii]->getCellModel())!=cm)
3042       throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::Merge1DGTUMeshes : all items must have the same geo type !");
3043   std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCoupling1DGTUMesh> > bb(sz);
3044   std::vector< const MEDCoupling1DGTUMesh * > aa(sz);
3045   int spaceDim=-3;
3046   for(std::size_t i=0;i<sz && spaceDim==-3;i++)
3047     {
3048       const MEDCoupling1DGTUMesh *cur=a[i];
3049       const DataArrayDouble *coo=cur->getCoords();
3050       if(coo)
3051         spaceDim=coo->getNumberOfComponents();
3052     }
3053   if(spaceDim==-3)
3054     throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::Merge1DGTUMeshes : no spaceDim specified ! unable to perform merge !");
3055   for(std::size_t i=0;i<sz;i++)
3056     {
3057       bb[i]=a[i]->buildSetInstanceFromThis(spaceDim);
3058       aa[i]=bb[i];
3059     }
3060   return Merge1DGTUMeshesLL(aa);
3061 }
3062
3063 /*!
3064  * \throw If presence of a null instance in the input vector \a a.
3065  * \throw If a is empty
3066  */
3067 MEDCoupling1DGTUMesh *MEDCoupling1DGTUMesh::Merge1DGTUMeshesOnSameCoords(std::vector<const MEDCoupling1DGTUMesh *>& a) throw(INTERP_KERNEL::Exception)
3068 {
3069   if(a.empty())
3070     throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::Merge1DGTUMeshesOnSameCoords : input array must be NON EMPTY !");
3071   std::vector<const MEDCoupling1DGTUMesh *>::const_iterator it=a.begin();
3072   if(!(*it))
3073     throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::Merge1DGTUMeshesOnSameCoords : null instance in the first element of input vector !");
3074   std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCoupling1DGTUMesh> > objs(a.size());
3075   std::vector<const DataArrayInt *> ncs(a.size()),ncis(a.size());
3076   int nbOfCells=(*it)->getNumberOfCells();
3077   const DataArrayDouble *coords=(*it)->getCoords();
3078   const INTERP_KERNEL::CellModel *cm=&((*it)->getCellModel());
3079   bool tmp;
3080   objs[0]=(*it)->copyWithNodalConnectivityPacked(tmp);
3081   ncs[0]=objs[0]->getNodalConnectivity(); ncis[0]=objs[0]->getNodalConnectivityIndex();
3082   it++;
3083   for(int i=1;it!=a.end();i++,it++)
3084     {
3085       if(!(*it))
3086         throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::Merge1DGTUMeshesOnSameCoords : presence of null instance !");
3087       if(cm!=&((*it)->getCellModel()))
3088         throw INTERP_KERNEL::Exception("Geometric types mismatches, Merge1DGTUMeshes impossible !");
3089       (*it)->getNumberOfCells();//to check that all is OK
3090       objs[i]=(*it)->copyWithNodalConnectivityPacked(tmp);
3091       ncs[i]=objs[i]->getNodalConnectivity(); ncis[i]=objs[i]->getNodalConnectivityIndex();
3092       if(coords!=(*it)->getCoords())
3093         throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::Merge1DGTUMeshesOnSameCoords : not lying on same coords !");
3094     }
3095   MEDCouplingAutoRefCountObjectPtr<MEDCoupling1DGTUMesh> ret(new MEDCoupling1DGTUMesh("merge",*cm));
3096   ret->setCoords(coords);
3097   ret->_conn=DataArrayInt::Aggregate(ncs);
3098   ret->_conn_indx=DataArrayInt::AggregateIndexes(ncis);
3099   return ret.retn();
3100 }
3101
3102 /*!
3103  * Assume that all instances in \a a are non null. If null it leads to a crash. That's why this method is assigned to be low level (LL)
3104  */
3105 MEDCoupling1DGTUMesh *MEDCoupling1DGTUMesh::Merge1DGTUMeshesLL(std::vector<const MEDCoupling1DGTUMesh *>& a) throw(INTERP_KERNEL::Exception)
3106 {
3107   if(a.empty())
3108     throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::Merge1DGTUMeshes : input array must be NON EMPTY !");
3109   std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCoupling1DGTUMesh> > objs(a.size());
3110   std::vector<const DataArrayInt *> ncs(a.size()),ncis(a.size());
3111   std::vector<const MEDCoupling1DGTUMesh *>::const_iterator it=a.begin();
3112   std::vector<int> nbNodesPerElt(a.size());
3113   int nbOfCells=(*it)->getNumberOfCells();
3114   bool tmp;
3115   objs[0]=(*it)->copyWithNodalConnectivityPacked(tmp);
3116   ncs[0]=objs[0]->getNodalConnectivity(); ncis[0]=objs[0]->getNodalConnectivityIndex();
3117   nbNodesPerElt[0]=0;
3118   int prevNbOfNodes=(*it)->getNumberOfNodes();
3119   const INTERP_KERNEL::CellModel *cm=&((*it)->getCellModel());
3120   it++;
3121   for(int i=1;it!=a.end();i++,it++)
3122     {
3123       if(cm!=&((*it)->getCellModel()))
3124         throw INTERP_KERNEL::Exception("Geometric types mismatches, Merge1DGTUMeshes impossible !");
3125       objs[i]=(*it)->copyWithNodalConnectivityPacked(tmp);
3126       ncs[i]=objs[i]->getNodalConnectivity(); ncis[i]=objs[i]->getNodalConnectivityIndex();
3127       nbOfCells+=(*it)->getNumberOfCells();
3128       nbNodesPerElt[i]=nbNodesPerElt[i-1]+prevNbOfNodes;
3129       prevNbOfNodes=(*it)->getNumberOfNodes();
3130     }
3131   std::vector<const MEDCouplingPointSet *> aps(a.size());
3132   std::copy(a.begin(),a.end(),aps.begin());
3133   MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> pts=MergeNodesArray(aps);
3134   MEDCouplingAutoRefCountObjectPtr<MEDCoupling1DGTUMesh> ret(new MEDCoupling1DGTUMesh("merge",*cm));
3135   ret->setCoords(pts);
3136   ret->_conn=AggregateNodalConnAndShiftNodeIds(ncs,nbNodesPerElt);
3137   ret->_conn_indx=DataArrayInt::AggregateIndexes(ncis);
3138   return ret.retn();
3139 }
3140
3141 MEDCoupling1DGTUMesh *MEDCoupling1DGTUMesh::buildSetInstanceFromThis(int spaceDim) const throw(INTERP_KERNEL::Exception)
3142 {
3143   MEDCouplingAutoRefCountObjectPtr<MEDCoupling1DGTUMesh> ret(new MEDCoupling1DGTUMesh(getName().c_str(),*_cm));
3144   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> tmp1,tmp2;
3145   const DataArrayInt *nodalConn(_conn),*nodalConnI(_conn_indx);
3146   if(!nodalConn)
3147     {
3148       tmp1=DataArrayInt::New(); tmp1->alloc(0,1);
3149     }
3150   else
3151     tmp1=_conn;
3152   ret->_conn=tmp1;
3153   //
3154   if(!nodalConnI)
3155     {
3156       tmp2=DataArrayInt::New(); tmp2->alloc(1,1); tmp2->setIJ(0,0,0);
3157     }
3158   else
3159     tmp2=_conn_indx;
3160   ret->_conn_indx=tmp2;
3161   //
3162   if(!_coords)
3163     {
3164       MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> coords=DataArrayDouble::New(); coords->alloc(0,spaceDim);
3165       ret->setCoords(coords);
3166     }
3167   else
3168     ret->setCoords(_coords);
3169   return ret.retn();
3170 }
3171
3172 /*!
3173  * This method aggregate the bbox of each cell and put it into bbox parameter.
3174  *
3175  * \return DataArrayDouble * - newly created object (to be managed by the caller) \a this number of cells tuples and 2*spacedim components.
3176  *
3177  * \throw If \a this is not fully set (coordinates and connectivity).
3178  * \throw If a cell in \a this has no valid nodeId.
3179  */
3180 DataArrayDouble *MEDCoupling1DGTUMesh::getBoundingBoxForBBTree() const
3181 {
3182   checkFullyDefined();
3183   int spaceDim(getSpaceDimension()),nbOfCells(getNumberOfCells()),nbOfNodes(getNumberOfNodes());
3184   MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> ret(DataArrayDouble::New()); ret->alloc(nbOfCells,2*spaceDim);
3185   double *bbox(ret->getPointer());
3186   for(int i=0;i<nbOfCells*spaceDim;i++)
3187     {
3188       bbox[2*i]=std::numeric_limits<double>::max();
3189       bbox[2*i+1]=-std::numeric_limits<double>::max();
3190     }
3191   const double *coordsPtr(_coords->getConstPointer());
3192   const int *conn(_conn->getConstPointer()),*connI(_conn_indx->getConstPointer());
3193   for(int i=0;i<nbOfCells;i++)
3194     {
3195       int offset=connI[i];
3196       int nbOfNodesForCell(connI[i+1]-offset),kk(0);
3197       for(int j=0;j<nbOfNodesForCell;j++)
3198         {
3199           int nodeId=conn[offset+j];
3200           if(nodeId>=0 && nodeId<nbOfNodes)
3201             {
3202               for(int k=0;k<spaceDim;k++)
3203                 {
3204                   bbox[2*spaceDim*i+2*k]=std::min(bbox[2*spaceDim*i+2*k],coordsPtr[spaceDim*nodeId+k]);
3205                   bbox[2*spaceDim*i+2*k+1]=std::max(bbox[2*spaceDim*i+2*k+1],coordsPtr[spaceDim*nodeId+k]);
3206                 }
3207               kk++;
3208             }
3209         }
3210       if(kk==0)
3211         {
3212           std::ostringstream oss; oss << "MEDCoupling1SGTUMesh::getBoundingBoxForBBTree : cell #" << i << " contains no valid nodeId !";
3213           throw INTERP_KERNEL::Exception(oss.str().c_str());
3214         }
3215     }
3216   return ret.retn();
3217 }
3218
3219 std::vector<int> MEDCoupling1DGTUMesh::BuildAPolygonFromParts(const std::vector< std::vector<int> >& parts) throw(INTERP_KERNEL::Exception)
3220 {
3221   std::vector<int> ret;
3222   if(parts.empty())
3223     return ret;
3224   ret.insert(ret.end(),parts[0].begin(),parts[0].end());
3225   int ref(ret.back());
3226   std::size_t sz(parts.size()),nbh(1);
3227   std::vector<bool> b(sz,true); b[0]=false;
3228   while(nbh<sz)
3229     {
3230       std::size_t i(0);
3231       for(;i<sz;i++) if(b[i] && parts[i].front()==ref) { ret.insert(ret.end(),parts[i].begin()+1,parts[i].end()); nbh++; break; }
3232       if(i<sz)
3233         ref=ret.back();
3234       else
3235         throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::BuildAPolygonFromParts : the input vector is not a part of a single polygon !");
3236     }
3237   if(ret.back()==ret.front())
3238     ret.pop_back();
3239   return ret;
3240 }
3241
3242 /*!
3243  * This method performs an aggregation of \a nodalConns (as DataArrayInt::Aggregate does) but in addition of that a shift is applied on the
3244  * values contained in \a nodalConns using corresponding offset specified in input \a offsetInNodeIdsPerElt.
3245  * But it also manage the values -1, that have a semantic in MEDCoupling1DGTUMesh class (separator for polyhedron).
3246  *
3247  * \param [in] nodalConns - a list of nodal connectivity arrays same size than \a offsetInNodeIdsPerElt.
3248  * \param [in] offsetInNodeIdsPerElt - a list of offsets to apply.
3249  * \return DataArrayInt * - A new object (to be managed by the caller) that is the result of the aggregation.
3250  * \throw If \a nodalConns or \a offsetInNodeIdsPerElt are empty.
3251  * \throw If \a nodalConns and \a offsetInNodeIdsPerElt have not the same size.
3252  * \throw If presence of null pointer in \a nodalConns.
3253  * \throw If presence of not allocated or array with not exactly one component in \a nodalConns.
3254  */
3255 DataArrayInt *MEDCoupling1DGTUMesh::AggregateNodalConnAndShiftNodeIds(const std::vector<const DataArrayInt *>& nodalConns, const std::vector<int>& offsetInNodeIdsPerElt) throw(INTERP_KERNEL::Exception)
3256 {
3257   std::size_t sz1(nodalConns.size()),sz2(offsetInNodeIdsPerElt.size());
3258   if(sz1!=sz2)
3259     throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::AggregateNodalConnAndShiftNodeIds : input vectors do not have the same size !");
3260   if(sz1==0)
3261     throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::AggregateNodalConnAndShiftNodeIds : empty vectors in input !");
3262   int nbOfTuples=0;
3263   for(std::vector<const DataArrayInt *>::const_iterator it=nodalConns.begin();it!=nodalConns.end();it++)
3264     {
3265       if(!(*it))
3266         throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::AggregateNodalConnAndShiftNodeIds : presence of null pointer in input vector !");
3267       if(!(*it)->isAllocated())
3268         throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::AggregateNodalConnAndShiftNodeIds : presence of non allocated array in input vector !");
3269       if((*it)->getNumberOfComponents()!=1)
3270         throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::AggregateNodalConnAndShiftNodeIds : presence of array with not exactly one component !");
3271       nbOfTuples+=(*it)->getNumberOfTuples();
3272     }
3273   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret=DataArrayInt::New(); ret->alloc(nbOfTuples,1);
3274   int *pt=ret->getPointer();
3275   int i=0;
3276   for(std::vector<const DataArrayInt *>::const_iterator it=nodalConns.begin();it!=nodalConns.end();it++,i++)
3277     {
3278       int curNbt=(*it)->getNumberOfTuples();
3279       const int *inPt=(*it)->begin();
3280       int offset=offsetInNodeIdsPerElt[i];
3281       for(int j=0;j<curNbt;j++,pt++)
3282         {
3283           if(inPt[j]!=-1)
3284             *pt=inPt[j]+offset;
3285           else
3286             *pt=-1;
3287         }
3288     }
3289   return ret.retn();
3290 }
3291
3292 MEDCoupling1DGTUMesh *MEDCoupling1DGTUMesh::New(const MEDCouplingUMesh *m) throw(INTERP_KERNEL::Exception)
3293 {
3294   if(!m)
3295     throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::New : input mesh is null !");
3296   std::set<INTERP_KERNEL::NormalizedCellType> gts(m->getAllGeoTypes());
3297   if(gts.size()!=1)
3298     throw INTERP_KERNEL::Exception("MEDCoupling1DGTUMesh::New : input mesh must have exactly one geometric type !");
3299   int geoType((int)*gts.begin());
3300   MEDCouplingAutoRefCountObjectPtr<MEDCoupling1DGTUMesh> ret(MEDCoupling1DGTUMesh::New(m->getName().c_str(),*gts.begin()));
3301   ret->setCoords(m->getCoords()); ret->setDescription(m->getDescription().c_str());
3302   int nbCells(m->getNumberOfCells());
3303   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> conn(DataArrayInt::New()),connI(DataArrayInt::New());
3304   conn->alloc(m->getMeshLength()-nbCells,1); connI->alloc(nbCells+1,1);
3305   int *c(conn->getPointer()),*ci(connI->getPointer()); *ci=0;
3306   const int *cin(m->getNodalConnectivity()->begin()),*ciin(m->getNodalConnectivityIndex()->begin());
3307   for(int i=0;i<nbCells;i++,ciin++,ci++)
3308     {
3309       if(cin[ciin[0]]==geoType)
3310         {
3311           if(ciin[1]-ciin[0]>=1)
3312             {
3313               c=std::copy(cin+ciin[0]+1,cin+ciin[1],c);
3314               ci[1]=ci[0]+ciin[1]-ciin[0]-1;
3315             }
3316           else
3317             {
3318               std::ostringstream oss; oss << "MEDCoupling1DGTUMesh::New(const MEDCouplingUMesh *m) : something is wrong in the input mesh at cell #" << i << " ! The size of cell is not >=0 !";
3319               throw INTERP_KERNEL::Exception(oss.str().c_str());
3320             }
3321         }
3322       else
3323         {
3324           std::ostringstream oss; oss << "MEDCoupling1DGTUMesh::New(const MEDCouplingUMesh *m) : something is wrong in the input mesh at cell #" << i << " ! The geometric type is not those expected !";
3325           throw INTERP_KERNEL::Exception(oss.str().c_str());
3326         }
3327     }
3328   ret->setNodalConnectivity(conn,connI);
3329   return ret.retn();
3330 }