Salome HOME
Unwarningization under Win.
[tools/medcoupling.git] / src / MEDCoupling / MEDCouplingUMesh.cxx
index fe6f8f71f00d178e2cd62e522d546a8a69d99d78..6bb2da093b7c8356683118aa3c83ca0635969d79 100644 (file)
@@ -33,7 +33,6 @@
 #include "InterpKernelMeshQuality.hxx"
 #include "InterpKernelCellSimplify.hxx"
 #include "InterpKernelGeo2DEdgeArcCircle.hxx"
-#include "MEDCouplingAutoRefCountObjectPtr.hxx"
 #include "InterpKernelAutoPtr.hxx"
 #include "InterpKernelGeo2DNode.hxx"
 #include "InterpKernelGeo2DEdgeLin.hxx"
@@ -89,6 +88,22 @@ MEDCouplingUMesh *MEDCouplingUMesh::clone(bool recDeepCpy) const
   return new MEDCouplingUMesh(*this,recDeepCpy);
 }
 
+/*!
+ * This method behaves mostly like MEDCouplingUMesh::deepCpy method, except that only nodal connectivity arrays are deeply copied.
+ * The coordinates are shared between \a this and the returned instance.
+ * 
+ * \return MEDCouplingUMesh * - A new object instance holding the copy of \a this (deep for connectivity, shallow for coordiantes)
+ * \sa MEDCouplingUMesh::deepCpy
+ */
+MEDCouplingPointSet *MEDCouplingUMesh::deepCpyConnectivityOnly() const throw(INTERP_KERNEL::Exception)
+{
+  checkConnectivityFullyDefined();
+  MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> ret=clone(false);
+  MEDCouplingAutoRefCountObjectPtr<DataArrayInt> c(getNodalConnectivity()->deepCpy()),ci(getNodalConnectivityIndex()->deepCpy());
+  ret->setConnectivity(c,ci);
+  return ret.retn();
+}
+
 void MEDCouplingUMesh::shallowCopyConnectivityFrom(const MEDCouplingPointSet *other) throw(INTERP_KERNEL::Exception)
 {
   if(!other)
@@ -100,14 +115,20 @@ void MEDCouplingUMesh::shallowCopyConnectivityFrom(const MEDCouplingPointSet *ot
   setConnectivity(otherC2->getNodalConnectivity(),otherC2->getNodalConnectivityIndex(),true);
 }
 
-std::size_t MEDCouplingUMesh::getHeapMemorySize() const
+std::size_t MEDCouplingUMesh::getHeapMemorySizeWithoutChildren() const
+{
+  std::size_t ret(MEDCouplingPointSet::getHeapMemorySizeWithoutChildren());
+  return ret;
+}
+
+std::vector<const BigMemoryObject *> MEDCouplingUMesh::getDirectChildren() const
 {
-  std::size_t ret=0;
+  std::vector<const BigMemoryObject *> ret(MEDCouplingPointSet::getDirectChildren());
   if(_nodal_connec)
-    ret+=_nodal_connec->getHeapMemorySize();
+    ret.push_back(_nodal_connec);
   if(_nodal_connec_index)
-    ret+=_nodal_connec_index->getHeapMemorySize();
-  return MEDCouplingPointSet::getHeapMemorySize()+ret;
+    ret.push_back(_nodal_connec_index);
+  return ret;
 }
 
 void MEDCouplingUMesh::updateTime() const
@@ -1375,6 +1396,7 @@ DataArrayInt *MEDCouplingUMesh::getNodeIdsInUse(int& nbrOfNodesInUse) const thro
  * So for pohyhedrons some nodes can be counted several times in the returned result.
  * 
  * \return a newly allocated array
+ * \sa MEDCouplingUMesh::computeEffectiveNbOfNodesPerCell
  */
 DataArrayInt *MEDCouplingUMesh::computeNbOfNodesPerCell() const throw(INTERP_KERNEL::Exception)
 {
@@ -1395,6 +1417,36 @@ DataArrayInt *MEDCouplingUMesh::computeNbOfNodesPerCell() const throw(INTERP_KER
   return ret.retn();
 }
 
+/*!
+ * This method computes effective number of nodes per cell. That is to say nodes appearing several times in nodal connectivity of a cell,
+ * will be counted only once here whereas it will be counted several times in MEDCouplingUMesh::computeNbOfNodesPerCell method.
+ *
+ * \return DataArrayInt * - new object to be deallocated by the caller.
+ * \sa MEDCouplingUMesh::computeNbOfNodesPerCell
+ */
+DataArrayInt *MEDCouplingUMesh::computeEffectiveNbOfNodesPerCell() const throw(INTERP_KERNEL::Exception)
+{
+  checkConnectivityFullyDefined();
+  int nbOfCells=getNumberOfCells();
+  MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret=DataArrayInt::New();
+  ret->alloc(nbOfCells,1);
+  int *retPtr=ret->getPointer();
+  const int *conn=getNodalConnectivity()->getConstPointer();
+  const int *connI=getNodalConnectivityIndex()->getConstPointer();
+  for(int i=0;i<nbOfCells;i++,retPtr++)
+    {
+      std::set<int> s(conn+connI[i]+1,conn+connI[i+1]);
+      if(conn[connI[i]]!=(int)INTERP_KERNEL::NORM_POLYHED)
+        *retPtr=(int)s.size();
+      else
+        {
+          s.erase(-1);
+          *retPtr=(int)s.size();
+        }
+    }
+  return ret.retn();
+}
+
 /*!
  * This method returns a newly allocated array containing this->getNumberOfCells() tuples and 1 component.
  * For each cell in \b this the number of faces constituting (entity of dimension this->getMeshDimension()-1) cell is computed.
@@ -1769,7 +1821,7 @@ bool MEDCouplingUMesh::areCellsIncludedIn(const MEDCouplingUMesh *other, int com
     }
   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> o2n=mesh->zipConnectivityTraducer(compType,nbOfCells);
   arr=o2n->substr(nbOfCells);
-  arr->setName(other->getName());
+  arr->setName(other->getName().c_str());
   int tmp;
   if(other->getNumberOfCells()==0)
     return true;
@@ -1813,7 +1865,7 @@ bool MEDCouplingUMesh::areCellsIncludedIn2(const MEDCouplingUMesh *other, DataAr
             }
         }
     }
-  arr2->setName(other->getName());
+  arr2->setName(other->getName().c_str());
   if(arr2->presenceOfValue(0))
     return false;
   arr=arr2.retn();
@@ -2783,7 +2835,7 @@ MEDCouplingUMesh *MEDCouplingUMesh::buildSetInstanceFromThis(int spaceDim) const
   int mdim=getMeshDimension();
   if(mdim<0)
     throw INTERP_KERNEL::Exception("MEDCouplingUMesh::buildSetInstanceFromThis : invalid mesh dimension ! Should be >= 0 !");
-  MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> ret=MEDCouplingUMesh::New(getName(),mdim);
+  MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> ret=MEDCouplingUMesh::New(getName().c_str(),mdim);
   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> tmp1,tmp2;
   bool needToCpyCT=true;
   if(!_nodal_connec)
@@ -3113,7 +3165,7 @@ MEDCouplingPointSet *MEDCouplingUMesh::buildPartOfMySelfKeepCoords2(int start, i
  */
 MEDCouplingPointSet *MEDCouplingUMesh::buildPartOfMySelfKeepCoords(const int *begin, const int *end) const
 {
-  checkFullyDefined();
+  checkConnectivityFullyDefined();
   int ncell=getNumberOfCells();
   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> ret=MEDCouplingUMesh::New();
   ret->_mesh_dim=_mesh_dim;
@@ -3830,13 +3882,13 @@ DataArrayDouble *MEDCouplingUMesh::distanceToPoints(const DataArrayDouble *pts,
   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret1=DataArrayInt::New(); ret1->alloc(nbOfPts,1);
   const int *nc=_nodal_connec->begin(),*ncI=_nodal_connec_index->begin(); const double *coords=_coords->begin();
   double *ret0Ptr=ret0->getPointer(); int *ret1Ptr=ret1->getPointer(); const double *ptsPtr=pts->begin();
-  std::vector<double> bbox;
-  getBoundingBoxForBBTree(bbox);
+  MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> bboxArr(getBoundingBoxForBBTree());
+  const double *bbox(bboxArr->begin());
   switch(spaceDim)
     {
     case 3:
       {
-        BBTreeDst<3> myTree(&bbox[0],0,0,nbCells);
+        BBTreeDst<3> myTree(bbox,0,0,nbCells);
         for(int i=0;i<nbOfPts;i++,ret0Ptr++,ret1Ptr++,ptsPtr+=3)
           {
             double x=std::numeric_limits<double>::max();
@@ -3849,7 +3901,7 @@ DataArrayDouble *MEDCouplingUMesh::distanceToPoints(const DataArrayDouble *pts,
       }
     case 2:
       {
-        BBTreeDst<2> myTree(&bbox[0],0,0,nbCells);
+        BBTreeDst<2> myTree(bbox,0,0,nbCells);
         for(int i=0;i<nbOfPts;i++,ret0Ptr++,ret1Ptr++,ptsPtr+=2)
           {
             double x=std::numeric_limits<double>::max();
@@ -3967,7 +4019,7 @@ int MEDCouplingUMesh::getCellContainingPoint(const double *pos, double eps) cons
  *          faster. 
  *  \param [in] pos - array of coordinates of the ball central point.
  *  \param [in] eps - ball radius.
- *  \param [in,out] elts - vector returning ids of the found cells. It is cleared
+ *  \param [out] elts - vector returning ids of the found cells. It is cleared
  *         before inserting ids.
  *  \throw If the coordinates array is not set.
  *  \throw If \a this->getMeshDimension() != \a this->getSpaceDimension().
@@ -3977,8 +4029,9 @@ int MEDCouplingUMesh::getCellContainingPoint(const double *pos, double eps) cons
  */
 void MEDCouplingUMesh::getCellsContainingPoint(const double *pos, double eps, std::vector<int>& elts) const
 {
-  std::vector<int> eltsIndex;
-  getCellsContainingPoints(pos,1,eps,elts,eltsIndex);
+  MEDCouplingAutoRefCountObjectPtr<DataArrayInt> eltsUg,eltsIndexUg;
+  getCellsContainingPoints(pos,1,eps,eltsUg,eltsIndexUg);
+  elts.clear(); elts.insert(elts.end(),eltsUg->begin(),eltsUg->end());
 }
 
 /// @cond INTERNAL
@@ -4107,13 +4160,12 @@ namespace ParaMEDMEM
 
 template<int SPACEDIM>
 void MEDCouplingUMesh::getCellsContainingPointsAlg(const double *coords, const double *pos, int nbOfPoints,
-                                                   double eps, std::vector<int>& elts, std::vector<int>& eltsIndex) const
+                                                   double eps, MEDCouplingAutoRefCountObjectPtr<DataArrayInt>& elts, MEDCouplingAutoRefCountObjectPtr<DataArrayInt>& eltsIndex) const
 {
-  std::vector<double> bbox;
-  eltsIndex.resize(nbOfPoints+1);
-  eltsIndex[0]=0;
-  elts.clear();
-  getBoundingBoxForBBTree(bbox);
+  elts=DataArrayInt::New(); eltsIndex=DataArrayInt::New(); eltsIndex->alloc(nbOfPoints+1,1); eltsIndex->setIJ(0,0,0); elts->alloc(0,1);
+  int *eltsIndexPtr(eltsIndex->getPointer());
+  MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> bboxArr(getBoundingBoxForBBTree());
+  const double *bbox(bboxArr->begin());
   int nbOfCells=getNumberOfCells();
   const int *conn=_nodal_connec->getConstPointer();
   const int *connI=_nodal_connec_index->getConstPointer();
@@ -4121,7 +4173,7 @@ void MEDCouplingUMesh::getCellsContainingPointsAlg(const double *coords, const d
   BBTree<SPACEDIM,int> myTree(&bbox[0],0,0,nbOfCells,-eps);
   for(int i=0;i<nbOfPoints;i++)
     {
-      eltsIndex[i+1]=eltsIndex[i];
+      eltsIndexPtr[i+1]=eltsIndexPtr[i];
       for(int j=0;j<SPACEDIM;j++)
         {
           bb[2*j]=pos[SPACEDIM*i+j];
@@ -4136,8 +4188,8 @@ void MEDCouplingUMesh::getCellsContainingPointsAlg(const double *coords, const d
                                                                                                (INTERP_KERNEL::NormalizedCellType)conn[connI[*iter]],
                                                                                                coords,conn+connI[*iter]+1,sz,eps))
             {
-              eltsIndex[i+1]++;
-              elts.push_back(*iter);
+              eltsIndexPtr[i+1]++;
+              elts->pushBackSilent(*iter);
             }
         }
     }
@@ -4151,8 +4203,8 @@ void MEDCouplingUMesh::getCellsContainingPointsAlg(const double *coords, const d
  *         this->getSpaceDimension() * \a nbOfPoints 
  *  \param [in] nbOfPoints - number of points to locate within \a this mesh.
  *  \param [in] eps - radius of balls (i.e. the precision).
- *  \param [in,out] elts - vector returning ids of found cells.
- *  \param [in,out] eltsIndex - an array, of length \a nbOfPoints + 1,
+ *  \param [out] elts - vector returning ids of found cells.
+ *  \param [out] eltsIndex - an array, of length \a nbOfPoints + 1,
  *         dividing cell ids in \a elts into groups each referring to one
  *         point. Its every element (except the last one) is an index pointing to the
  *         first id of a group of cells. For example cells in contact with the *i*-th
@@ -4168,7 +4220,7 @@ void MEDCouplingUMesh::getCellsContainingPointsAlg(const double *coords, const d
  *  \ref  py_mcumesh_getCellsContainingPoints "Here is a Python example".
  */
 void MEDCouplingUMesh::getCellsContainingPoints(const double *pos, int nbOfPoints, double eps,
-                                                std::vector<int>& elts, std::vector<int>& eltsIndex) const
+                                                MEDCouplingAutoRefCountObjectPtr<DataArrayInt>& elts, MEDCouplingAutoRefCountObjectPtr<DataArrayInt>& eltsIndex) const
 {
   int spaceDim=getSpaceDimension();
   int mDim=getMeshDimension();
@@ -5192,6 +5244,7 @@ void MEDCouplingUMesh::tessellate2DCurve(double eps) throw(INTERP_KERNEL::Except
  *          and \a this->getMeshDimension() != 3. 
  *  \throw If \a policy is not one of the four discussed above.
  *  \throw If the nodal connectivity of cells is not defined.
+ * \sa MEDCouplingUMesh::tetrahedrize
  */
 DataArrayInt *MEDCouplingUMesh::simplexize(int policy) throw(INTERP_KERNEL::Exception)
 {
@@ -6128,36 +6181,49 @@ MEDCouplingFieldDouble *MEDCouplingUMesh::getSkewField() const throw(INTERP_KERN
 
 /*!
  * This method aggregate the bbox of each cell and put it into bbox parameter.
- * \param bbox out parameter of size 2*spacedim*nbOfcells.
+ * 
+ * \return DataArrayDouble * - newly created object (to be managed by the caller) \a this number of cells tuples and 2*spacedim components.
+ * 
+ * \throw If \a this is not fully set (coordinates and connectivity).
+ * \throw If a cell in \a this has no valid nodeId.
  */
-void MEDCouplingUMesh::getBoundingBoxForBBTree(std::vector<double>& bbox) const
+DataArrayDouble *MEDCouplingUMesh::getBoundingBoxForBBTree() const
 {
-  int spaceDim=getSpaceDimension();
-  int nbOfCells=getNumberOfCells();
-  bbox.resize(2*nbOfCells*spaceDim);
+  checkFullyDefined();
+  int spaceDim(getSpaceDimension()),nbOfCells(getNumberOfCells()),nbOfNodes(getNumberOfNodes());
+  MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> ret(DataArrayDouble::New()); ret->alloc(nbOfCells,2*spaceDim);
+  double *bbox(ret->getPointer());
   for(int i=0;i<nbOfCells*spaceDim;i++)
     {
       bbox[2*i]=std::numeric_limits<double>::max();
       bbox[2*i+1]=-std::numeric_limits<double>::max();
     }
-  const double *coordsPtr=_coords->getConstPointer();
-  const int *conn=_nodal_connec->getConstPointer();
-  const int *connI=_nodal_connec_index->getConstPointer();
+  const double *coordsPtr(_coords->getConstPointer());
+  const int *conn(_nodal_connec->getConstPointer()),*connI(_nodal_connec_index->getConstPointer());
   for(int i=0;i<nbOfCells;i++)
     {
       int offset=connI[i]+1;
-      int nbOfNodesForCell=connI[i+1]-offset;
+      int nbOfNodesForCell(connI[i+1]-offset),kk(0);
       for(int j=0;j<nbOfNodesForCell;j++)
         {
           int nodeId=conn[offset+j];
-          if(nodeId>=0)
-            for(int k=0;k<spaceDim;k++)
-              {
-                bbox[2*spaceDim*i+2*k]=std::min(bbox[2*spaceDim*i+2*k],coordsPtr[spaceDim*nodeId+k]);
-                bbox[2*spaceDim*i+2*k+1]=std::max(bbox[2*spaceDim*i+2*k+1],coordsPtr[spaceDim*nodeId+k]);
-              }
+          if(nodeId>=0 && nodeId<nbOfNodes)
+            {
+              for(int k=0;k<spaceDim;k++)
+                {
+                  bbox[2*spaceDim*i+2*k]=std::min(bbox[2*spaceDim*i+2*k],coordsPtr[spaceDim*nodeId+k]);
+                  bbox[2*spaceDim*i+2*k+1]=std::max(bbox[2*spaceDim*i+2*k+1],coordsPtr[spaceDim*nodeId+k]);
+                }
+              kk++;
+            }
+        }
+      if(kk==0)
+        {
+          std::ostringstream oss; oss << "MEDCouplingUMesh::getBoundingBoxForBBTree : cell #" << i << " contains no valid nodeId !";
+          throw INTERP_KERNEL::Exception(oss.str().c_str());
         }
     }
+  return ret.retn();
 }
 
 /// @cond INTERNAL
@@ -6192,7 +6258,7 @@ namespace ParaMEDMEMImpl
  * This method returns in the same format as code (see MEDCouplingUMesh::checkTypeConsistencyAndContig or MEDCouplingUMesh::splitProfilePerType) how
  * \a this is composed in cell types.
  * The returned array is of size 3*n where n is the number of different types present in \a this. 
- * For every k in [0,n] ret[3*k+2]==0 because it has no sense here. 
+ * For every k in [0,n] ret[3*k+2]==-1 because it has no sense here. 
  * This parameter is kept only for compatibility with other methode listed above.
  */
 std::vector<int> MEDCouplingUMesh::getDistributionOfTypes() const throw(INTERP_KERNEL::Exception)
@@ -6203,7 +6269,7 @@ std::vector<int> MEDCouplingUMesh::getDistributionOfTypes() const throw(INTERP_K
   const int *work=connI;
   int nbOfCells=getNumberOfCells();
   std::size_t n=getAllTypes().size();
-  std::vector<int> ret(3*n,0); //ret[3*k+2]==0 because it has no sense here
+  std::vector<int> ret(3*n,-1); //ret[3*k+2]==-1 because it has no sense here
   std::set<INTERP_KERNEL::NormalizedCellType> types;
   for(std::size_t i=0;work!=connI+nbOfCells;i++)
     {
@@ -6655,7 +6721,7 @@ DataArrayInt *MEDCouplingUMesh::rearrange2ConsecutiveCellTypes()
  */
 std::vector<MEDCouplingUMesh *> MEDCouplingUMesh::splitByType() const
 {
-  checkFullyDefined();
+  checkConnectivityFullyDefined();
   const int *conn=_nodal_connec->getConstPointer();
   const int *connI=_nodal_connec_index->getConstPointer();
   int nbOfCells=getNumberOfCells();
@@ -6693,7 +6759,7 @@ MEDCoupling1GTUMesh *MEDCouplingUMesh::convertIntoSingleGeoTypeMesh() const thro
     if(_types.size()!=1)
     throw INTERP_KERNEL::Exception("MEDCouplingUMesh::convertIntoSingleGeoTypeMesh : current mesh does not contain exactly one geometric type !");
   INTERP_KERNEL::NormalizedCellType typ=*_types.begin();
-  MEDCouplingAutoRefCountObjectPtr<MEDCoupling1GTUMesh> ret=MEDCoupling1GTUMesh::New(getName(),typ);
+  MEDCouplingAutoRefCountObjectPtr<MEDCoupling1GTUMesh> ret=MEDCoupling1GTUMesh::New(getName().c_str(),typ);
   ret->setCoords(getCoords());
   MEDCoupling1SGTUMesh *retC=dynamic_cast<MEDCoupling1SGTUMesh *>((MEDCoupling1GTUMesh*)ret);
   if(retC)
@@ -7392,7 +7458,7 @@ MEDCouplingUMesh *MEDCouplingUMesh::FuseUMeshesOnSameCoords(const std::vector<co
       tmp->alloc(curNbOfCells,1);
       std::copy(o2nPtr+offset,o2nPtr+offset+curNbOfCells,tmp->getPointer());
       offset+=curNbOfCells;
-      tmp->setName(meshes[i]->getName());
+      tmp->setName(meshes[i]->getName().c_str());
       corr[i]=tmp;
     }
   return ret.retn();
@@ -8009,7 +8075,7 @@ void MEDCouplingUMesh::FillInCompact3DMode(int spaceDim, int nbOfNodesInCell, co
     throw INTERP_KERNEL::Exception("MEDCouplingUMesh::FillInCompact3DMode : Invalid spaceDim specified : must be 2 or 3 !");
 }
 
-void MEDCouplingUMesh::writeVTKLL(std::ostream& ofs, const std::string& cellData, const std::string& pointData) const throw(INTERP_KERNEL::Exception)
+void MEDCouplingUMesh::writeVTKLL(std::ostream& ofs, const std::string& cellData, const std::string& pointData, DataArrayByte *byteData) const throw(INTERP_KERNEL::Exception)
 {
   int nbOfCells=getNumberOfCells();
   if(nbOfCells<=0)
@@ -8023,11 +8089,11 @@ void MEDCouplingUMesh::writeVTKLL(std::ostream& ofs, const std::string& cellData
   ofs << "      </CellData>\n";
   ofs << "      <Points>\n";
   if(getSpaceDimension()==3)
-    _coords->writeVTK(ofs,8,"Points");
+    _coords->writeVTK(ofs,8,"Points",byteData);
   else
     {
       MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> coo=_coords->changeNbOfComponents(3,0.);
-      coo->writeVTK(ofs,8,"Points");
+      coo->writeVTK(ofs,8,"Points",byteData);
     }
   ofs << "      </Points>\n";
   ofs << "      <Cells>\n";
@@ -8058,12 +8124,12 @@ void MEDCouplingUMesh::writeVTKLL(std::ostream& ofs, const std::string& cellData
         }
     }
   types->transformWithIndArr(PARAMEDMEM2VTKTYPETRADUCER,PARAMEDMEM2VTKTYPETRADUCER+INTERP_KERNEL::NORM_MAXTYPE);
-  types->writeVTK(ofs,8,"UInt8","types");
-  offsets->writeVTK(ofs,8,"Int32","offsets");
+  types->writeVTK(ofs,8,"UInt8","types",byteData);
+  offsets->writeVTK(ofs,8,"Int32","offsets",byteData);
   if(szFaceOffsets!=0)
     {//presence of Polyhedra
       connectivity->reAlloc(szConn);
-      faceoffsets->writeVTK(ofs,8,"Int32","faceoffsets");
+      faceoffsets->writeVTK(ofs,8,"Int32","faceoffsets",byteData);
       MEDCouplingAutoRefCountObjectPtr<DataArrayInt> faces=DataArrayInt::New(); faces->alloc(szFaceOffsets,1);
       w1=faces->getPointer();
       for(int i=0;i<nbOfCells;i++)
@@ -8080,9 +8146,9 @@ void MEDCouplingUMesh::writeVTKLL(std::ostream& ofs, const std::string& cellData
                 w6=w5+1;
               }
           }
-      faces->writeVTK(ofs,8,"Int32","faces");
+      faces->writeVTK(ofs,8,"Int32","faces",byteData);
     }
-  connectivity->writeVTK(ofs,8,"Int32","connectivity");
+  connectivity->writeVTK(ofs,8,"Int32","connectivity",byteData);
   ofs << "      </Cells>\n";
   ofs << "    </Piece>\n";
   ofs << "  </" << getVTKDataSetType() << ">\n";
@@ -8193,7 +8259,6 @@ void MEDCouplingUMesh::BuildIntersecting2DCellsFromEdges(double eps, const MEDCo
                                                          std::vector<double>& addCoordsQuadratic, std::vector<int>& cr, std::vector<int>& crI, std::vector<int>& cNb1, std::vector<int>& cNb2)
 {
   static const int SPACEDIM=2;
-  std::vector<double> bbox1,bbox2;
   const double *coo1=m1->getCoords()->getConstPointer();
   const int *conn1=m1->getNodalConnectivity()->getConstPointer();
   const int *connI1=m1->getNodalConnectivityIndex()->getConstPointer();
@@ -8203,15 +8268,15 @@ void MEDCouplingUMesh::BuildIntersecting2DCellsFromEdges(double eps, const MEDCo
   const int *connI2=m2->getNodalConnectivityIndex()->getConstPointer();
   int offset2=offset1+m2->getNumberOfNodes();
   int offset3=offset2+((int)addCoords.size())/2;
-  m1->getBoundingBoxForBBTree(bbox1);
-  m2->getBoundingBoxForBBTree(bbox2);
-  BBTree<SPACEDIM,int> myTree(&bbox2[0],0,0,m2->getNumberOfCells(),eps);
+  MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> bbox1Arr(m1->getBoundingBoxForBBTree()),bbox2Arr(m2->getBoundingBoxForBBTree());
+  const double *bbox1(bbox1Arr->begin()),*bbox2(bbox2Arr->begin());
+  BBTree<SPACEDIM,int> myTree(bbox2,0,0,m2->getNumberOfCells(),eps);
   int ncell1=m1->getNumberOfCells();
   crI.push_back(0);
   for(int i=0;i<ncell1;i++)
     {
       std::vector<int> candidates2;
-      myTree.getIntersectingElems(&bbox1[i*2*SPACEDIM],candidates2);
+      myTree.getIntersectingElems(bbox1+i*2*SPACEDIM,candidates2);
       std::map<INTERP_KERNEL::Node *,int> mapp;
       std::map<int,INTERP_KERNEL::Node *> mappRev;
       INTERP_KERNEL::QuadraticPolygon pol1;
@@ -8286,22 +8351,21 @@ void MEDCouplingUMesh::IntersectDescending2DMeshes(const MEDCouplingUMesh *m1, c
   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> dd9(m1Desc),dd10(m2Desc);
   const int *c1=m1Desc->getNodalConnectivity()->getConstPointer();
   const int *ci1=m1Desc->getNodalConnectivityIndex()->getConstPointer();
-  std::vector<double> bbox1,bbox2;
-  m1Desc->getBoundingBoxForBBTree(bbox1);
-  m2Desc->getBoundingBoxForBBTree(bbox2);
+  MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> bbox1Arr(m1Desc->getBoundingBoxForBBTree()),bbox2Arr(m2Desc->getBoundingBoxForBBTree());
+  const double *bbox1(bbox1Arr->begin()),*bbox2(bbox2Arr->begin());
   int ncell1=m1Desc->getNumberOfCells();
   int ncell2=m2Desc->getNumberOfCells();
   intersectEdge1.resize(ncell1);
   colinear2.resize(ncell2);
   subDiv2.resize(ncell2);
-  BBTree<SPACEDIM,int> myTree(&bbox2[0],0,0,m2Desc->getNumberOfCells(),-eps);
+  BBTree<SPACEDIM,int> myTree(bbox2,0,0,m2Desc->getNumberOfCells(),-eps);
   std::vector<int> candidates1(1);
   int offset1=m1->getNumberOfNodes();
   int offset2=offset1+m2->getNumberOfNodes();
   for(int i=0;i<ncell1;i++)
     {
       std::vector<int> candidates2;
-      myTree.getIntersectingElems(&bbox1[i*2*SPACEDIM],candidates2);
+      myTree.getIntersectingElems(bbox1+i*2*SPACEDIM,candidates2);
       if(!candidates2.empty())
         {
           std::map<INTERP_KERNEL::Node *,int> map1,map2;
@@ -8787,7 +8851,7 @@ void MEDCouplingUMesh::ExtractFromIndexedArrays2(int idsOfSelectStart, int idsOf
   int *work=arrIo->getPointer();
   *work++=0;
   int lgth=0;
-  for(std::size_t i=0;i<sz;i++,work++,idsIt+=idsOfSelectStep)
+  for(int i=0;i<sz;i++,work++,idsIt+=idsOfSelectStep)
     {
       if(idsIt>=0 && idsIt<nbOfGrps)
         lgth+=arrIndxPtr[idsIt+1]-arrIndxPtr[idsIt];
@@ -8808,7 +8872,7 @@ void MEDCouplingUMesh::ExtractFromIndexedArrays2(int idsOfSelectStart, int idsOf
   arro->alloc(lgth,1);
   work=arro->getPointer();
   idsIt=idsOfSelectStart;
-  for(std::size_t i=0;i<sz;i++,idsIt+=idsOfSelectStep)
+  for(int i=0;i<sz;i++,idsIt+=idsOfSelectStep)
     {
       if(arrIndxPtr[idsIt]>=0 && arrIndxPtr[idsIt+1]<=maxSizeOfArr)
         work=std::copy(arrInPtr+arrIndxPtr[idsIt],arrInPtr+arrIndxPtr[idsIt+1],work);
@@ -9159,7 +9223,7 @@ MEDCouplingUMesh *MEDCouplingUMesh::buildSpreadZonesWithPoly() const throw(INTER
   std::vector<DataArrayInt *> partition=partitionBySpreadZone();
   std::vector< MEDCouplingAutoRefCountObjectPtr<DataArrayInt> > partitionAuto; partitionAuto.reserve(partition.size());
   std::copy(partition.begin(),partition.end(),std::back_insert_iterator<std::vector< MEDCouplingAutoRefCountObjectPtr<DataArrayInt> > >(partitionAuto));
-  MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> ret=MEDCouplingUMesh::New(getName(),mdim);
+  MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> ret=MEDCouplingUMesh::New(getName().c_str(),mdim);
   ret->setCoords(getCoords());
   ret->allocateCells((int)partition.size());
   //
@@ -9193,7 +9257,6 @@ MEDCouplingUMesh *MEDCouplingUMesh::buildSpreadZonesWithPoly() const throw(INTER
  */
 std::vector<DataArrayInt *> MEDCouplingUMesh::partitionBySpreadZone() const throw(INTERP_KERNEL::Exception)
 {
-  //#if 0
   int nbOfCellsCur=getNumberOfCells();
   std::vector<DataArrayInt *> ret;
   if(nbOfCellsCur<=0)
@@ -9213,36 +9276,6 @@ std::vector<DataArrayInt *> MEDCouplingUMesh::partitionBySpreadZone() const thro
   for(std::vector< MEDCouplingAutoRefCountObjectPtr<DataArrayInt> >::iterator it=ret2.begin();it!=ret2.end();it++)
     ret.push_back((*it).retn());
   return ret;
-  //#endif
-#if 0
-  int nbOfCellsCur=getNumberOfCells();
-  DataArrayInt *neigh=0,*neighI=0;
-  computeNeighborsOfCells(neigh,neighI);
-  MEDCouplingAutoRefCountObjectPtr<DataArrayInt> neighAuto(neigh),neighIAuto(neighI);
-  MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ids=DataArrayInt::New(); ids->alloc(nbOfCellsCur,1); ids->iota();
-  std::vector<DataArrayInt *> ret;
-  std::vector< MEDCouplingAutoRefCountObjectPtr<DataArrayInt> > ret2;
-  while(nbOfCellsCur>0)
-    {
-      MEDCouplingAutoRefCountObjectPtr<DataArrayInt> tmp=MEDCouplingUMesh::ComputeSpreadZoneGradually(neighAuto,neighIAuto);
-      MEDCouplingAutoRefCountObjectPtr<DataArrayInt> tmp3=tmp->buildComplement(nbOfCellsCur);
-      MEDCouplingAutoRefCountObjectPtr<DataArrayInt> tmp2=ids->selectByTupleId(tmp->begin(),tmp->end());
-      ret2.push_back(tmp2);  ret.push_back(tmp2);
-      nbOfCellsCur=tmp3->getNumberOfTuples();
-      if(nbOfCellsCur>0)
-        {
-          ids=ids->selectByTupleId(tmp3->begin(),tmp3->end());
-          MEDCouplingUMesh::ExtractFromIndexedArrays(tmp3->begin(),tmp3->end(),neighAuto,neighIAuto,neigh,neighI);
-          neighAuto=neigh;
-          neighIAuto=neighI;
-          MEDCouplingAutoRefCountObjectPtr<DataArrayInt> renum=tmp3->invertArrayN2O2O2N(nbOfCellsCur+tmp->getNumberOfTuples());
-          neighAuto->transformWithIndArr(renum->begin(),renum->end());
-        }
-    }
-  for(std::vector<DataArrayInt *>::const_iterator it=ret.begin();it!=ret.end();it++)
-    (*it)->incrRef();
-  return ret;
-#endif
 }
 
 /*!
@@ -9269,6 +9302,78 @@ DataArrayInt *MEDCouplingUMesh::ComputeRangesFromTypeDistribution(const std::vec
   return ret.retn();
 }
 
+/*!
+ * This method expects that \a this a 3D mesh (spaceDim=3 and meshDim=3) with all coordinates and connectivities set.
+ * All cells in \a this are expected to be linear 3D cells.
+ * This method will split **all** 3D cells in \a this into INTERP_KERNEL::NORM_TETRA4 cells and put them in the returned mesh.
+ * It leads to an increase to number of cells.
+ * This method contrary to MEDCouplingUMesh::simplexize can append coordinates in \a this to perform its work.
+ * The \a nbOfAdditionalPoints returned value informs about it. If > 0, the coordinates array in returned mesh will have \a nbOfAdditionalPoints 
+ * more tuples (nodes) than in \a this. Anyway, all the nodes in \a this (with the same order) will be in the returned mesh.
+ *
+ * \param [in] policy - the policy of splitting that must be in (PLANAR_FACE_5, PLANAR_FACE_6, GENERAL_24, GENERAL_48). The policy will be used only for INTERP_KERNEL::NORM_HEXA8 cells.
+ *                      For all other cells, the splitting policy will be ignored.
+ * \param [out] nbOfAdditionalPoints - number of nodes added to \c this->_coords. If > 0 a new coordinates object will be constructed result of the aggregation of the old one and the new points added. 
+ * \param [out] n2oCells - A new instance of DataArrayInt holding, for each new cell,
+ *          an id of old cell producing it. The caller is to delete this array using
+ *         decrRef() as it is no more needed.
+ * \return MEDCoupling1SGTUMesh * - the mesh containing only INTERP_KERNEL::NORM_TETRA4 cells.
+ * 
+ * \throw If \a this is not a 3D mesh (spaceDim==3 and meshDim==3).
+ * \throw If \a this is not fully constituted with linear 3D cells.
+ * \sa MEDCouplingUMesh::simplexize
+ */
+MEDCoupling1SGTUMesh *MEDCouplingUMesh::tetrahedrize(int policy, DataArrayInt *& n2oCells, int& nbOfAdditionalPoints) const throw(INTERP_KERNEL::Exception)
+{
+  INTERP_KERNEL::SplittingPolicy pol((INTERP_KERNEL::SplittingPolicy)policy);
+  checkConnectivityFullyDefined();
+  if(getMeshDimension()!=3 || getSpaceDimension()!=3)
+    throw INTERP_KERNEL::Exception("MEDCouplingUMesh::tetrahedrize : only available for mesh with meshdim == 3 and spacedim == 3 !");
+  int nbOfCells(getNumberOfCells()),nbNodes(getNumberOfNodes());
+  MEDCouplingAutoRefCountObjectPtr<MEDCoupling1SGTUMesh> ret0(MEDCoupling1SGTUMesh::New(getName().c_str(),INTERP_KERNEL::NORM_TETRA4));
+  MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret(DataArrayInt::New()); ret->alloc(nbOfCells,1);
+  int *retPt(ret->getPointer());
+  MEDCouplingAutoRefCountObjectPtr<DataArrayInt> newConn(DataArrayInt::New()); newConn->alloc(0,1);
+  MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> addPts(DataArrayDouble::New()); addPts->alloc(0,1);
+  const int *oldc(_nodal_connec->begin());
+  const int *oldci(_nodal_connec_index->begin());
+  const double *coords(_coords->begin());
+  for(int i=0;i<nbOfCells;i++,oldci++,retPt++)
+    {
+      std::vector<int> a; std::vector<double> b;
+      INTERP_KERNEL::SplitIntoTetras(pol,(INTERP_KERNEL::NormalizedCellType)oldc[oldci[0]],oldc+oldci[0]+1,oldc+oldci[1],coords,a,b);
+      std::size_t nbOfTet(a.size()/4); *retPt=(int)nbOfTet;
+      const int *aa(&a[0]);
+      if(!b.empty())
+        {
+          for(std::vector<int>::iterator it=a.begin();it!=a.end();it++)
+            if(*it<0)
+              *it=(-(*(it))-1+nbNodes);
+          addPts->insertAtTheEnd(b.begin(),b.end());
+          nbNodes+=(int)b.size()/3;
+        }
+      for(std::size_t j=0;j<nbOfTet;j++,aa+=4)
+        newConn->insertAtTheEnd(aa,aa+4);
+    }
+  if(!addPts->empty())
+    {
+      addPts->rearrange(3);
+      nbOfAdditionalPoints=addPts->getNumberOfTuples();
+      addPts=DataArrayDouble::Aggregate(getCoords(),addPts);
+      ret0->setCoords(addPts);
+    }
+  else
+    {
+      nbOfAdditionalPoints=0;
+      ret0->setCoords(getCoords());
+    }
+  ret0->setNodalConnectivity(newConn);
+  //
+  ret->computeOffsets2();
+  n2oCells=ret->buildExplicitArrOfSliceOnScaledArr(0,nbOfCells,1);
+  return ret0.retn();
+}
+
 MEDCouplingUMeshCellIterator::MEDCouplingUMeshCellIterator(MEDCouplingUMesh *mesh):_mesh(mesh),_cell(new MEDCouplingUMeshCell(mesh)),
                                                                                    _own_cell(true),_cell_id(-1),_nb_cell(0)
 {