Salome HOME
cherry-pick between a669 and 19a96 of V7_main
[tools/medcoupling.git] / src / MEDCoupling / MEDCouplingUMesh.cxx
index 26fd3a45b54d8edf1c81afb8007ed6a37d9ed824..10cbc6bf23fba7c299d7b14269552cdab7be349e 100644 (file)
 #include "BBTreeDst.txx"
 #include "SplitterTetra.hxx"
 #include "DirectedBoundingBox.hxx"
+#include "InterpKernelMatrixTools.hxx"
 #include "InterpKernelMeshQuality.hxx"
 #include "InterpKernelCellSimplify.hxx"
 #include "InterpKernelGeo2DEdgeArcCircle.hxx"
-#include "MEDCouplingAutoRefCountObjectPtr.hxx"
 #include "InterpKernelAutoPtr.hxx"
 #include "InterpKernelGeo2DNode.hxx"
 #include "InterpKernelGeo2DEdgeLin.hxx"
@@ -58,7 +58,7 @@ MEDCouplingUMesh *MEDCouplingUMesh::New()
   return new MEDCouplingUMesh;
 }
 
-MEDCouplingUMesh *MEDCouplingUMesh::New(const char *meshName, int meshDim)
+MEDCouplingUMesh *MEDCouplingUMesh::New(const std::string& meshName, int meshDim)
 {
   MEDCouplingUMesh *ret=new MEDCouplingUMesh;
   ret->setName(meshName);
@@ -89,7 +89,23 @@ MEDCouplingUMesh *MEDCouplingUMesh::clone(bool recDeepCpy) const
   return new MEDCouplingUMesh(*this,recDeepCpy);
 }
 
-void MEDCouplingUMesh::shallowCopyConnectivityFrom(const MEDCouplingPointSet *other) throw(INTERP_KERNEL::Exception)
+/*!
+ * 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
+{
+  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)
 {
   if(!other)
     throw INTERP_KERNEL::Exception("MEDCouplingUMesh::shallowCopyConnectivityFrom : input pointer is null !");
@@ -100,14 +116,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
@@ -140,7 +162,7 @@ MEDCouplingUMesh::MEDCouplingUMesh():_mesh_dim(-2),_nodal_connec(0),_nodal_conne
  *  \throw If the connectivity index data array has more than one component.
  *  \throw If the connectivity index data array has a named component.
  */
-void MEDCouplingUMesh::checkCoherency() const throw(INTERP_KERNEL::Exception)
+void MEDCouplingUMesh::checkCoherency() const
 {
   if(_mesh_dim<-1)
    throw INTERP_KERNEL::Exception("No mesh dimension specified !");
@@ -193,7 +215,7 @@ void MEDCouplingUMesh::checkCoherency() const throw(INTERP_KERNEL::Exception)
  *  \throw If number of nodes defining an element does not correspond to the type of element.
  *  \throw If the nodal connectivity includes an invalid node id.
  */
-void MEDCouplingUMesh::checkCoherency1(double eps) const throw(INTERP_KERNEL::Exception)
+void MEDCouplingUMesh::checkCoherency1(double eps) const
 {
   checkCoherency();
   if(_mesh_dim==-1)
@@ -265,7 +287,7 @@ void MEDCouplingUMesh::checkCoherency1(double eps) const throw(INTERP_KERNEL::Ex
  *  \throw If number of nodes defining an element does not correspond to the type of element.
  *  \throw If the nodal connectivity includes an invalid node id.
  */
-void MEDCouplingUMesh::checkCoherency2(double eps) const throw(INTERP_KERNEL::Exception)
+void MEDCouplingUMesh::checkCoherency2(double eps) const
 {
   checkCoherency1(eps);
 }
@@ -292,8 +314,10 @@ void MEDCouplingUMesh::setMeshDimension(int meshDim)
  *
  *  \param [in] nbOfCells - estimation of the number of cell \a this mesh will contain.
  *
+ *  \if ENABLE_EXAMPLES
  *  \ref medcouplingcppexamplesUmeshStdBuild1 "Here is a C++ example".<br>
  *  \ref medcouplingpyexamplesUmeshStdBuild1 "Here is a Python example".
+ *  \endif
  */
 void MEDCouplingUMesh::allocateCells(int nbOfCells)
 {
@@ -323,10 +347,12 @@ void MEDCouplingUMesh::allocateCells(int nbOfCells)
  *  \param [in] size - number of nodes constituting this cell.
  *  \param [in] nodalConnOfCell - the connectivity of the cell to add.
  * 
+ *  \if ENABLE_EXAMPLES
  *  \ref medcouplingcppexamplesUmeshStdBuild1 "Here is a C++ example".<br>
  *  \ref medcouplingpyexamplesUmeshStdBuild1 "Here is a Python example".
+ *  \endif
  */
-void MEDCouplingUMesh::insertNextCell(INTERP_KERNEL::NormalizedCellType type, int size, const int *nodalConnOfCell) throw(INTERP_KERNEL::Exception)
+void MEDCouplingUMesh::insertNextCell(INTERP_KERNEL::NormalizedCellType type, int size, const int *nodalConnOfCell)
 {
   const INTERP_KERNEL::CellModel& cm=INTERP_KERNEL::CellModel::GetCellModel(type);
   if(_nodal_connec_index==0)
@@ -359,8 +385,10 @@ void MEDCouplingUMesh::insertNextCell(INTERP_KERNEL::NormalizedCellType type, in
  * Compacts data arrays to release unused memory. This method is to be called after
  * finishing cell insertion using \a this->insertNextCell().
  * 
+ *  \if ENABLE_EXAMPLES
  *  \ref medcouplingcppexamplesUmeshStdBuild1 "Here is a C++ example".<br>
  *  \ref medcouplingpyexamplesUmeshStdBuild1 "Here is a Python example".
+ *  \endif
  */
 void MEDCouplingUMesh::finishInsertingCells()
 {
@@ -386,7 +414,7 @@ MEDCouplingUMeshCellIterator *MEDCouplingUMesh::cellIterator()
  * In this case MEDCouplingUMesh::sortCellsInMEDFileFrmt or MEDCouplingUMesh::rearrange2ConsecutiveCellTypes methods for example can be called before invoking this method.
  * Useful for python users.
  */
-MEDCouplingUMeshCellByTypeEntry *MEDCouplingUMesh::cellsByType() throw(INTERP_KERNEL::Exception)
+MEDCouplingUMeshCellByTypeEntry *MEDCouplingUMesh::cellsByType()
 {
   if(!checkConsecutiveCellTypes())
     throw INTERP_KERNEL::Exception("MEDCouplingUMesh::cellsByType : this mesh is not sorted by type !");
@@ -397,17 +425,44 @@ MEDCouplingUMeshCellByTypeEntry *MEDCouplingUMesh::cellsByType() throw(INTERP_KE
  * Returns a set of all cell types available in \a this mesh.
  * \return std::set<INTERP_KERNEL::NormalizedCellType> - the set of cell types.
  * \warning this method does not throw any exception even if \a this is not defined.
+ * \sa MEDCouplingUMesh::getAllGeoTypesSorted
  */
 std::set<INTERP_KERNEL::NormalizedCellType> MEDCouplingUMesh::getAllGeoTypes() const
 {
   return _types;
 }
 
+/*!
+ * This method returns the sorted list of geometric types in \a this.
+ * Sorted means in the same order than the cells in \a this. A single entry in return vector means the maximal chunk of consecutive cells in \a this
+ * having the same geometric type. So a same geometric type can appear more than once if the cells are not sorted per geometric type.
+ *
+ * \throw if connectivity in \a this is not correctly defined.
+ *  
+ * \sa MEDCouplingMesh::getAllGeoTypes
+ */
+std::vector<INTERP_KERNEL::NormalizedCellType> MEDCouplingUMesh::getAllGeoTypesSorted() const
+{
+  std::vector<INTERP_KERNEL::NormalizedCellType> ret;
+  checkConnectivityFullyDefined();
+  int nbOfCells(getNumberOfCells());
+  if(nbOfCells==0)
+    return ret;
+  if(getMeshLength()<1)
+    throw INTERP_KERNEL::Exception("MEDCouplingUMesh::getAllGeoTypesSorted : the connectivity in this seems invalid !");
+  const int *c(_nodal_connec->begin()),*ci(_nodal_connec_index->begin());
+  ret.push_back((INTERP_KERNEL::NormalizedCellType)c[*ci++]);
+  for(int i=1;i<nbOfCells;i++,ci++)
+    if(ret.back()!=((INTERP_KERNEL::NormalizedCellType)c[*ci]))
+      ret.push_back((INTERP_KERNEL::NormalizedCellType)c[*ci]);
+  return ret;
+}
+
 /*!
  * This method is a method that compares \a this and \a other.
  * This method compares \b all attributes, even names and component names.
  */
-bool MEDCouplingUMesh::isEqualIfNotWhy(const MEDCouplingMesh *other, double prec, std::string& reason) const throw(INTERP_KERNEL::Exception)
+bool MEDCouplingUMesh::isEqualIfNotWhy(const MEDCouplingMesh *other, double prec, std::string& reason) const
 {
   if(!other)
     throw INTERP_KERNEL::Exception("MEDCouplingUMesh::isEqualIfNotWhy : input other pointer is null !");
@@ -508,7 +563,7 @@ bool MEDCouplingUMesh::isEqualWithoutConsideringStr(const MEDCouplingMesh *other
  *  \param [in] prec - the precision used to compare nodes of the two meshes.
  *  \throw If the two meshes do not match.
  */
-void MEDCouplingUMesh::checkFastEquivalWith(const MEDCouplingMesh *other, double prec) const throw(INTERP_KERNEL::Exception)
+void MEDCouplingUMesh::checkFastEquivalWith(const MEDCouplingMesh *other, double prec) const
 {
  MEDCouplingPointSet::checkFastEquivalWith(other,prec);
  const MEDCouplingUMesh *otherC=dynamic_cast<const MEDCouplingUMesh *>(other);
@@ -534,10 +589,12 @@ void MEDCouplingUMesh::checkFastEquivalWith(const MEDCouplingMesh *other, double
  * \throw If the coordinates array is not set.
  * \throw If the nodal connectivity of cells is not defined.
  * 
+ * \if ENABLE_EXAMPLES
  * \ref cpp_mcumesh_getReverseNodalConnectivity "Here is a C++ example".<br>
  * \ref  py_mcumesh_getReverseNodalConnectivity "Here is a Python example".
+ * \endif
  */
-void MEDCouplingUMesh::getReverseNodalConnectivity(DataArrayInt *revNodal, DataArrayInt *revNodalIndx) const throw(INTERP_KERNEL::Exception)
+void MEDCouplingUMesh::getReverseNodalConnectivity(DataArrayInt *revNodal, DataArrayInt *revNodalIndx) const
 {
   checkFullyDefined();
   int nbOfNodes=getNumberOfNodes();
@@ -666,11 +723,13 @@ private:
  *  \throw If \a desc == NULL || \a descIndx == NULL || \a revDesc == NULL || \a
  *         revDescIndx == NULL.
  * 
+ *  \if ENABLE_EXAMPLES
  *  \ref cpp_mcumesh_buildDescendingConnectivity "Here is a C++ example".<br>
  *  \ref  py_mcumesh_buildDescendingConnectivity "Here is a Python example".
+ *  \endif
  * \sa buildDescendingConnectivity2()
  */
-MEDCouplingUMesh *MEDCouplingUMesh::buildDescendingConnectivity(DataArrayInt *desc, DataArrayInt *descIndx, DataArrayInt *revDesc, DataArrayInt *revDescIndx) const throw(INTERP_KERNEL::Exception)
+MEDCouplingUMesh *MEDCouplingUMesh::buildDescendingConnectivity(DataArrayInt *desc, DataArrayInt *descIndx, DataArrayInt *revDesc, DataArrayInt *revDescIndx) const
 {
   return buildDescendingConnectivityGen<MinusOneSonsGenerator>(desc,descIndx,revDesc,revDescIndx,MEDCouplingFastNbrer);
 }
@@ -682,7 +741,7 @@ MEDCouplingUMesh *MEDCouplingUMesh::buildDescendingConnectivity(DataArrayInt *de
  * This method returns 4 arrays and a mesh as MEDCouplingUMesh::buildDescendingConnectivity does.
  * \sa MEDCouplingUMesh::buildDescendingConnectivity
  */
-MEDCouplingUMesh *MEDCouplingUMesh::explode3DMeshTo1D(DataArrayInt *desc, DataArrayInt *descIndx, DataArrayInt *revDesc, DataArrayInt *revDescIndx) const throw(INTERP_KERNEL::Exception)
+MEDCouplingUMesh *MEDCouplingUMesh::explode3DMeshTo1D(DataArrayInt *desc, DataArrayInt *descIndx, DataArrayInt *revDesc, DataArrayInt *revDescIndx) const
 {
   checkFullyDefined();
   if(getMeshDimension()!=3)
@@ -735,11 +794,13 @@ MEDCouplingUMesh *MEDCouplingUMesh::explode3DMeshTo1D(DataArrayInt *desc, DataAr
  *  \throw If \a desc == NULL || \a descIndx == NULL || \a revDesc == NULL || \a
  *         revDescIndx == NULL.
  * 
+ *  \if ENABLE_EXAMPLES
  *  \ref cpp_mcumesh_buildDescendingConnectivity2 "Here is a C++ example".<br>
  *  \ref  py_mcumesh_buildDescendingConnectivity2 "Here is a Python example".
+ *  \endif
  * \sa buildDescendingConnectivity()
  */
-MEDCouplingUMesh *MEDCouplingUMesh::buildDescendingConnectivity2(DataArrayInt *desc, DataArrayInt *descIndx, DataArrayInt *revDesc, DataArrayInt *revDescIndx) const throw(INTERP_KERNEL::Exception)
+MEDCouplingUMesh *MEDCouplingUMesh::buildDescendingConnectivity2(DataArrayInt *desc, DataArrayInt *descIndx, DataArrayInt *revDesc, DataArrayInt *revDescIndx) const
 {
   return buildDescendingConnectivityGen<MinusOneSonsGenerator>(desc,descIndx,revDesc,revDescIndx,MEDCouplingOrientationSensitiveNbrer);
 }
@@ -754,7 +815,7 @@ MEDCouplingUMesh *MEDCouplingUMesh::buildDescendingConnectivity2(DataArrayInt *d
  *                        parameter allows to select the right part in this array. The number of tuples is equal to the last values in \b neighborsIndx.
  * \param [out] neighborsIndx is an array of size this->getNumberOfCells()+1 newly allocated and should be dealt by the caller. This arrays allow to use the first output parameter \b neighbors.
  */
-void MEDCouplingUMesh::computeNeighborsOfCells(DataArrayInt *&neighbors, DataArrayInt *&neighborsIndx) const throw(INTERP_KERNEL::Exception)
+void MEDCouplingUMesh::computeNeighborsOfCells(DataArrayInt *&neighbors, DataArrayInt *&neighborsIndx) const
 {
   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> desc=DataArrayInt::New();
   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> descIndx=DataArrayInt::New();
@@ -818,7 +879,7 @@ void MEDCouplingUMesh::ComputeNeighborsOfCellsAdv(const DataArrayInt *desc, cons
  * For speed reasons no check of this will be done.
  */
 template<class SonsGenerator>
-MEDCouplingUMesh *MEDCouplingUMesh::buildDescendingConnectivityGen(DataArrayInt *desc, DataArrayInt *descIndx, DataArrayInt *revDesc, DataArrayInt *revDescIndx, DimM1DescNbrer nbrer) const throw(INTERP_KERNEL::Exception)
+MEDCouplingUMesh *MEDCouplingUMesh::buildDescendingConnectivityGen(DataArrayInt *desc, DataArrayInt *descIndx, DataArrayInt *revDesc, DataArrayInt *revDescIndx, DimM1DescNbrer nbrer) const
 {
   if(!desc || !descIndx || !revDesc || !revDescIndx)
     throw INTERP_KERNEL::Exception("MEDCouplingUMesh::buildDescendingConnectivityGen : present of a null pointer in input !");
@@ -830,7 +891,7 @@ MEDCouplingUMesh *MEDCouplingUMesh::buildDescendingConnectivityGen(DataArrayInt
   const int *conn=_nodal_connec->getConstPointer();
   const int *connIndex=_nodal_connec_index->getConstPointer();
   std::string name="Mesh constituent of "; name+=getName();
-  MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> ret=MEDCouplingUMesh::New(name.c_str(),getMeshDimension()-SonsGenerator::DELTA);
+  MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> ret=MEDCouplingUMesh::New(name,getMeshDimension()-SonsGenerator::DELTA);
   ret->setCoords(getCoords());
   ret->allocateCells(2*nbOfCells);
   descIndx->alloc(nbOfCells+1,1);
@@ -957,8 +1018,10 @@ struct MEDCouplingAccVisit
  *  \throw If the nodal connectivity of cells is node defined.
  *  \throw If dimension of \a this mesh is not either 2 or 3.
  *
+ *  \if ENABLE_EXAMPLES
  *  \ref cpp_mcumesh_convertToPolyTypes "Here is a C++ example".<br>
  *  \ref  py_mcumesh_convertToPolyTypes "Here is a Python example".
+ *  \endif
  */
 void MEDCouplingUMesh::convertToPolyTypes(const int *cellIdsToConvertBg, const int *cellIdsToConvertEnd)
 {
@@ -966,7 +1029,7 @@ void MEDCouplingUMesh::convertToPolyTypes(const int *cellIdsToConvertBg, const i
   int dim=getMeshDimension();
   if(dim<2 || dim>3)
     throw INTERP_KERNEL::Exception("Invalid mesh dimension : must be 2 or 3 !");
-  int nbOfCells=getNumberOfCells();
+  int nbOfCells(getNumberOfCells());
   if(dim==2)
     {
       const int *connIndex=_nodal_connec_index->getConstPointer();
@@ -991,47 +1054,50 @@ void MEDCouplingUMesh::convertToPolyTypes(const int *cellIdsToConvertBg, const i
     }
   else
     {
-      int *connIndex=_nodal_connec_index->getPointer();
-      int connIndexLgth=_nodal_connec_index->getNbOfElems();
-      const int *connOld=_nodal_connec->getConstPointer();
-      int connOldLgth=_nodal_connec->getNbOfElems();
-      std::vector<int> connNew(connOld,connOld+connOldLgth);
+      int *connIndex(_nodal_connec_index->getPointer());
+      const int *connOld(_nodal_connec->getConstPointer());
+      MEDCouplingAutoRefCountObjectPtr<DataArrayInt> connNew(DataArrayInt::New()),connNewI(DataArrayInt::New()); connNew->alloc(0,1); connNewI->alloc(1,1); connNewI->setIJ(0,0,0);
+      std::vector<bool> toBeDone(nbOfCells,false);
       for(const int *iter=cellIdsToConvertBg;iter!=cellIdsToConvertEnd;iter++)
         {
           if(*iter>=0 && *iter<nbOfCells)
+            toBeDone[*iter]=true;
+          else
             {
-              int pos=connIndex[*iter];
-              int posP1=connIndex[(*iter)+1];
-              int lgthOld=posP1-pos-1;
-              const INTERP_KERNEL::CellModel& cm=INTERP_KERNEL::CellModel::GetCellModel((INTERP_KERNEL::NormalizedCellType)connNew[pos]);
-              connNew[pos]=INTERP_KERNEL::NORM_POLYHED;
-              unsigned nbOfFaces=cm.getNumberOfSons2(&connNew[pos+1],lgthOld);
-              int *tmp=new int[nbOfFaces*lgthOld];
-              int *work=tmp;
-              for(int j=0;j<(int)nbOfFaces;j++)
+              std::ostringstream oss; oss << "MEDCouplingUMesh::convertToPolyTypes : On rank #" << std::distance(cellIdsToConvertBg,iter) << " value is " << *iter << " which is not";
+              oss << " in range [0," << nbOfCells << ") !";
+              throw INTERP_KERNEL::Exception(oss.str().c_str());
+            }
+        }
+      for(int cellId=0;cellId<nbOfCells;cellId++)
+        {
+          int pos(connIndex[cellId]),posP1(connIndex[cellId+1]);
+          int lgthOld(posP1-pos-1);
+          if(toBeDone[cellId])
+            {
+              const INTERP_KERNEL::CellModel& cm=INTERP_KERNEL::CellModel::GetCellModel((INTERP_KERNEL::NormalizedCellType)connOld[pos]);
+              unsigned nbOfFaces(cm.getNumberOfSons2(connOld+pos+1,lgthOld));
+              int *tmp(new int[nbOfFaces*lgthOld+1]);
+              int *work=tmp; *work++=INTERP_KERNEL::NORM_POLYHED;
+              for(unsigned j=0;j<nbOfFaces;j++)
                 {
                   INTERP_KERNEL::NormalizedCellType type;
-                  unsigned offset=cm.fillSonCellNodalConnectivity2(j,&connNew[pos+1],lgthOld,work,type);
+                  unsigned offset=cm.fillSonCellNodalConnectivity2(j,connOld+pos+1,lgthOld,work,type);
                   work+=offset;
                   *work++=-1;
                 }
-              std::size_t newLgth=std::distance(tmp,work)-1;
-              std::size_t delta=newLgth-lgthOld;
-              std::transform(connIndex+(*iter)+1,connIndex+connIndexLgth,connIndex+(*iter)+1,std::bind2nd(std::plus<int>(),delta));
-              connNew.insert(connNew.begin()+posP1,tmp+lgthOld,tmp+newLgth);
-              std::copy(tmp,tmp+lgthOld,connNew.begin()+pos+1);
+              std::size_t newLgth(std::distance(tmp,work)-1);//-1 for last -1
+              connNew->pushBackValsSilent(tmp,tmp+newLgth);
+              connNewI->pushBackSilent(connNewI->back()+(int)newLgth);
               delete [] tmp;
             }
           else
             {
-              std::ostringstream oss; oss << "MEDCouplingUMesh::convertToPolyTypes : On rank #" << std::distance(cellIdsToConvertBg,iter) << " value is " << *iter << " which is not";
-              oss << " in range [0," << nbOfCells << ") !";
-              throw INTERP_KERNEL::Exception(oss.str().c_str());
+              connNew->pushBackValsSilent(connOld+pos,connOld+posP1);
+              connNewI->pushBackSilent(connNewI->back()+posP1-pos);
             }
         }
-      _nodal_connec->alloc((int)connNew.size(),1);
-      int *newConnPtr=_nodal_connec->getPointer();
-      std::copy(connNew.begin(),connNew.end(),newConnPtr);
+      setConnectivity(connNew,connNewI,false);//false because computeTypes called just behind.
     }
   computeTypes();
 }
@@ -1077,10 +1143,12 @@ void MEDCouplingUMesh::convertAllToPoly()
  *  \throw If \a this mesh contains polyhedrons with the valid connectivity.
  *  \throw If \a this mesh contains polyhedrons with odd number of nodes.
  *
+ *  \if ENABLE_EXAMPLES
  *  \ref cpp_mcumesh_arePolyhedronsNotCorrectlyOriented "Here is a C++ example".<br>
  *  \ref  py_mcumesh_arePolyhedronsNotCorrectlyOriented "Here is a Python example".
+ *  \endif
  */
-void MEDCouplingUMesh::convertExtrudedPolyhedra() throw(INTERP_KERNEL::Exception)
+void MEDCouplingUMesh::convertExtrudedPolyhedra()
 {
   checkFullyDefined();
   if(getMeshDimension()!=3 || getSpaceDimension()!=3)
@@ -1234,7 +1302,7 @@ bool MEDCouplingUMesh::unPolyze()
  * \param [in] eps is a relative precision that allows to establish if some 3D plane are coplanar or not. This epsilon is used to recenter around origin to have maximal 
  *             precision.
  */
-void MEDCouplingUMesh::simplifyPolyhedra(double eps) throw(INTERP_KERNEL::Exception)
+void MEDCouplingUMesh::simplifyPolyhedra(double eps)
 {
   checkFullyDefined();
   if(getMeshDimension()!=3 || getSpaceDimension()!=3)
@@ -1273,7 +1341,7 @@ void MEDCouplingUMesh::simplifyPolyhedra(double eps) throw(INTERP_KERNEL::Except
  * \return a newly allocated DataArrayInt sorted ascendingly of fetched node ids.
  * \sa MEDCouplingUMesh::getNodeIdsInUse
  */
-DataArrayInt *MEDCouplingUMesh::computeFetchedNodeIds() const throw(INTERP_KERNEL::Exception)
+DataArrayInt *MEDCouplingUMesh::computeFetchedNodeIds() const
 {
   checkConnectivityFullyDefined();
   int nbOfCells=getNumberOfCells();
@@ -1303,7 +1371,7 @@ DataArrayInt *MEDCouplingUMesh::computeFetchedNodeIds() const throw(INTERP_KERNE
  * \param [in,out] nodeIdsInUse an array of size typically equal to nbOfNodes.
  * \sa MEDCouplingUMesh::getNodeIdsInUse
  */
-void MEDCouplingUMesh::computeNodeIdsAlg(std::vector<bool>& nodeIdsInUse) const throw(INTERP_KERNEL::Exception)
+void MEDCouplingUMesh::computeNodeIdsAlg(std::vector<bool>& nodeIdsInUse) const
 {
   int nbOfNodes=(int)nodeIdsInUse.size();
   int nbOfCells=getNumberOfCells();
@@ -1336,11 +1404,13 @@ void MEDCouplingUMesh::computeNodeIdsAlg(std::vector<bool>& nodeIdsInUse) const
  *  \throw If the nodal connectivity of cells is not defined.
  *  \throw If the nodal connectivity includes an invalid id.
  *
+ *  \if ENABLE_EXAMPLES
  *  \ref cpp_mcumesh_getNodeIdsInUse "Here is a C++ example".<br>
  *  \ref  py_mcumesh_getNodeIdsInUse "Here is a Python example".
+ *  \endif
  * \sa computeNodeIdsAlg()
  */
-DataArrayInt *MEDCouplingUMesh::getNodeIdsInUse(int& nbrOfNodesInUse) const throw(INTERP_KERNEL::Exception)
+DataArrayInt *MEDCouplingUMesh::getNodeIdsInUse(int& nbrOfNodesInUse) const
 {
   nbrOfNodesInUse=-1;
   int nbOfNodes=getNumberOfNodes();
@@ -1375,8 +1445,9 @@ 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)
+DataArrayInt *MEDCouplingUMesh::computeNbOfNodesPerCell() const
 {
   checkConnectivityFullyDefined();
   int nbOfCells=getNumberOfCells();
@@ -1395,13 +1466,43 @@ 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
+{
+  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.
  * 
  * \return a newly allocated array
  */
-DataArrayInt *MEDCouplingUMesh::computeNbOfFacesPerCell() const throw(INTERP_KERNEL::Exception)
+DataArrayInt *MEDCouplingUMesh::computeNbOfFacesPerCell() const
 {
   checkConnectivityFullyDefined();
   int nbOfCells=getNumberOfCells();
@@ -1429,10 +1530,12 @@ DataArrayInt *MEDCouplingUMesh::computeNbOfFacesPerCell() const throw(INTERP_KER
  *  \throw If the nodal connectivity of cells is not defined.
  *  \throw If the nodal connectivity includes an invalid id.
  *
+ *  \if ENABLE_EXAMPLES
  *  \ref cpp_mcumesh_zipCoordsTraducer "Here is a C++ example".<br>
  *  \ref  py_mcumesh_zipCoordsTraducer "Here is a Python example".
+ *  \endif
  */
-DataArrayInt *MEDCouplingUMesh::zipCoordsTraducer() throw(INTERP_KERNEL::Exception)
+DataArrayInt *MEDCouplingUMesh::zipCoordsTraducer()
 {
   return MEDCouplingPointSet::zipCoordsTraducer();
 }
@@ -1645,7 +1748,7 @@ bool MEDCouplingUMesh::AreCellsEqualInPool(const std::vector<int>& candidates, i
  * \return the correspondance array old to new in a newly allocated array.
  * 
  */
-void MEDCouplingUMesh::findCommonCells(int compType, int startCellId, DataArrayInt *& commonCellsArr, DataArrayInt *& commonCellsIArr) const throw(INTERP_KERNEL::Exception)
+void MEDCouplingUMesh::findCommonCells(int compType, int startCellId, DataArrayInt *& commonCellsArr, DataArrayInt *& commonCellsIArr) const
 {
   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> revNodal=DataArrayInt::New(),revNodalI=DataArrayInt::New();
   getReverseNodalConnectivity(revNodal,revNodalI);
@@ -1750,12 +1853,14 @@ void MEDCouplingUMesh::FindCommonCellsAlg(int compType, int startCellId, const D
  *  \return bool - \c true if all cells of \a other mesh are present in the \a this
  *         mesh.
  *
+ *  \if ENABLE_EXAMPLES
  *  \ref cpp_mcumesh_areCellsIncludedIn "Here is a C++ example".<br>
  *  \ref  py_mcumesh_areCellsIncludedIn "Here is a Python example".
+ *  \endif
  *  \sa checkDeepEquivalOnSameNodesWith()
  *  \sa checkGeoEquivalWith()
  */
-bool MEDCouplingUMesh::areCellsIncludedIn(const MEDCouplingUMesh *other, int compType, DataArrayInt *& arr) const throw(INTERP_KERNEL::Exception)
+bool MEDCouplingUMesh::areCellsIncludedIn(const MEDCouplingUMesh *other, int compType, DataArrayInt *& arr) const
 {
   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> mesh=MergeUMeshesOnSameCoords(this,other);
   int nbOfCells=getNumberOfCells();
@@ -1785,7 +1890,7 @@ bool MEDCouplingUMesh::areCellsIncludedIn(const MEDCouplingUMesh *other, int com
  * \param arr is an output parameter that returns a \b newly created instance. This array is of size 'other->getNumberOfCells()'.
  * \return If \a other is fully included in 'this 'true is returned. If not false is returned.
  */
-bool MEDCouplingUMesh::areCellsIncludedIn2(const MEDCouplingUMesh *other, DataArrayInt *& arr) const throw(INTERP_KERNEL::Exception)
+bool MEDCouplingUMesh::areCellsIncludedIn2(const MEDCouplingUMesh *other, DataArrayInt *& arr) const
 {
   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> mesh=MergeUMeshesOnSameCoords(this,other);
   DataArrayInt *commonCells=0,*commonCellsI=0;
@@ -1844,7 +1949,7 @@ MEDCouplingPointSet *MEDCouplingUMesh::mergeMyselfWithOnSameCoords(const MEDCoup
  * \warning This method modifies can generate an unstructured mesh whose cells are not sorted by geometric type order.
  * In view of the MED file writing, a renumbering of cells of returned unstructured mesh (using MEDCouplingUMesh::sortCellsInMEDFileFrmt) should be necessary.
  */
-MEDCouplingPointSet *MEDCouplingUMesh::buildPartOfMySelf2(int start, int end, int step, bool keepCoords) const throw(INTERP_KERNEL::Exception)
+MEDCouplingPointSet *MEDCouplingUMesh::buildPartOfMySelf2(int start, int end, int step, bool keepCoords) const
 {
   if(getMeshDimension()!=-1)
     return MEDCouplingPointSet::buildPartOfMySelf2(start,end,step,keepCoords);
@@ -1878,8 +1983,10 @@ MEDCouplingPointSet *MEDCouplingUMesh::buildPartOfMySelf2(int start, int end, in
  *  \throw If the nodal connectivity of cells is not defined.
  *  \throw If any cell id in the array \a begin is not valid.
  *
+ *  \if ENABLE_EXAMPLES
  *  \ref cpp_mcumesh_buildPartOfMySelf "Here is a C++ example".<br>
  *  \ref  py_mcumesh_buildPartOfMySelf "Here is a Python example".
+ *  \endif
  */
 MEDCouplingPointSet *MEDCouplingUMesh::buildPartOfMySelf(const int *begin, const int *end, bool keepCoords) const
 {
@@ -1908,7 +2015,7 @@ MEDCouplingPointSet *MEDCouplingUMesh::buildPartOfMySelf(const int *begin, const
  * \param [in] otherOnSameCoordsThanThis an another mesh with same meshdimension than \b this with exactly the same number of cells than cell ids list in [\b cellIdsBg, \b cellIdsEnd ).
  *             Coordinate pointer of \b this and those of \b otherOnSameCoordsThanThis must be the same
  */
-void MEDCouplingUMesh::setPartOfMySelf(const int *cellIdsBg, const int *cellIdsEnd, const MEDCouplingUMesh& otherOnSameCoordsThanThis) throw(INTERP_KERNEL::Exception)
+void MEDCouplingUMesh::setPartOfMySelf(const int *cellIdsBg, const int *cellIdsEnd, const MEDCouplingUMesh& otherOnSameCoordsThanThis)
 {
   checkConnectivityFullyDefined();
   otherOnSameCoordsThanThis.checkConnectivityFullyDefined();
@@ -1957,7 +2064,7 @@ void MEDCouplingUMesh::setPartOfMySelf(const int *cellIdsBg, const int *cellIdsE
     }
 }
 
-void MEDCouplingUMesh::setPartOfMySelf2(int start, int end, int step, const MEDCouplingUMesh& otherOnSameCoordsThanThis) throw(INTERP_KERNEL::Exception)
+void MEDCouplingUMesh::setPartOfMySelf2(int start, int end, int step, const MEDCouplingUMesh& otherOnSameCoordsThanThis)
 {
   checkConnectivityFullyDefined();
   otherOnSameCoordsThanThis.checkConnectivityFullyDefined();
@@ -2066,8 +2173,10 @@ void MEDCouplingUMesh::fillCellIdsToKeepFromNodeIds(const int *begin, const int
  *  \throw If the nodal connectivity of cells is not defined.
  *  \throw If any node id in \a begin is not valid.
  *
+ *  \if ENABLE_EXAMPLES
  *  \ref cpp_mcumesh_buildFacePartOfMySelfNode "Here is a C++ example".<br>
  *  \ref  py_mcumesh_buildFacePartOfMySelfNode "Here is a Python example".
+ *  \endif
  */
 MEDCouplingPointSet *MEDCouplingUMesh::buildFacePartOfMySelfNode(const int *begin, const int *end, bool fullyIn) const
 {
@@ -2089,8 +2198,10 @@ MEDCouplingPointSet *MEDCouplingUMesh::buildFacePartOfMySelfNode(const int *begi
  *  \throw If the coordinates array is not set.
  *  \throw If the nodal connectivity of cells is not defined.
  *
+ *  \if ENABLE_EXAMPLES
  *  \ref cpp_mcumesh_buildBoundaryMesh "Here is a C++ example".<br>
  *  \ref  py_mcumesh_buildBoundaryMesh "Here is a Python example".
+ *  \endif
  */
 MEDCouplingPointSet *MEDCouplingUMesh::buildBoundaryMesh(bool keepCoords) const
 {
@@ -2119,7 +2230,7 @@ MEDCouplingPointSet *MEDCouplingUMesh::buildBoundaryMesh(bool keepCoords) const
  * A cell is detected to be on boundary if it contains one or more than one face having only one father.
  * This method makes the assumption that \a this is fully defined (coords,connectivity). If not an exception will be thrown. 
  */
-DataArrayInt *MEDCouplingUMesh::findCellIdsOnBoundary() const throw(INTERP_KERNEL::Exception)
+DataArrayInt *MEDCouplingUMesh::findCellIdsOnBoundary() const
 {
   checkFullyDefined();
   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> desc=DataArrayInt::New();
@@ -2171,7 +2282,7 @@ DataArrayInt *MEDCouplingUMesh::findCellIdsOnBoundary() const throw(INTERP_KERNE
  * \param [out] cellIdsRk1 a newly allocated array containing cells ids of s1+s2 \b into \b cellIdsRk0 subset. To get absolute ids of s1+s2 simply invoke
  *              cellIdsRk1->transformWithIndArr(cellIdsRk0->begin(),cellIdsRk0->end());
  */
-void MEDCouplingUMesh::findCellIdsLyingOn(const MEDCouplingUMesh& otherDimM1OnSameCoords, DataArrayInt *&cellIdsRk0, DataArrayInt *&cellIdsRk1) const throw(INTERP_KERNEL::Exception)
+void MEDCouplingUMesh::findCellIdsLyingOn(const MEDCouplingUMesh& otherDimM1OnSameCoords, DataArrayInt *&cellIdsRk0, DataArrayInt *&cellIdsRk1) const
 {
   if(getCoords()!=otherDimM1OnSameCoords.getCoords())
     throw INTERP_KERNEL::Exception("MEDCouplingUMesh::findCellIdsLyingOn : coordinates pointer are not the same ! Use tryToShareSameCoords method !");
@@ -2219,7 +2330,7 @@ void MEDCouplingUMesh::findCellIdsLyingOn(const MEDCouplingUMesh& otherDimM1OnSa
  * 
  * \return a newly allocated mesh lying on the same coordinates than \b this. The caller has to deal with returned mesh.
  */
-MEDCouplingUMesh *MEDCouplingUMesh::computeSkin() const throw(INTERP_KERNEL::Exception)
+MEDCouplingUMesh *MEDCouplingUMesh::computeSkin() const
 {
   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> desc=DataArrayInt::New();
   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> descIndx=DataArrayInt::New();
@@ -2241,8 +2352,10 @@ MEDCouplingUMesh *MEDCouplingUMesh::computeSkin() const throw(INTERP_KERNEL::Exc
  *  \throw If the coordinates array is not set.
  *  \throw If the nodal connectivity of cells is node defined.
  *
+ *  \if ENABLE_EXAMPLES
  *  \ref cpp_mcumesh_findBoundaryNodes "Here is a C++ example".<br>
  *  \ref  py_mcumesh_findBoundaryNodes "Here is a Python example".
+ *  \endif
  */
 DataArrayInt *MEDCouplingUMesh::findBoundaryNodes() const
 {
@@ -2250,7 +2363,7 @@ DataArrayInt *MEDCouplingUMesh::findBoundaryNodes() const
   return skin->computeFetchedNodeIds();
 }
 
-MEDCouplingUMesh *MEDCouplingUMesh::buildUnstructured() const throw(INTERP_KERNEL::Exception)
+MEDCouplingUMesh *MEDCouplingUMesh::buildUnstructured() const
 {
   incrRef();
   return const_cast<MEDCouplingUMesh *>(this);
@@ -2327,7 +2440,7 @@ void MEDCouplingUMesh::findNodesToDuplicate(const MEDCouplingUMesh& otherDimM1On
  * \param [in] nodeIdsToDuplicateBg begin of node ids (included) to be duplicated in connectivity only
  * \param [in] nodeIdsToDuplicateEnd end of node ids (excluded) to be duplicated in connectivity only
  */
-void MEDCouplingUMesh::duplicateNodes(const int *nodeIdsToDuplicateBg, const int *nodeIdsToDuplicateEnd) throw(INTERP_KERNEL::Exception)
+void MEDCouplingUMesh::duplicateNodes(const int *nodeIdsToDuplicateBg, const int *nodeIdsToDuplicateEnd)
 {
   int nbOfNodes=getNumberOfNodes();
   duplicateNodesInCoords(nodeIdsToDuplicateBg,nodeIdsToDuplicateEnd);
@@ -2344,8 +2457,10 @@ void MEDCouplingUMesh::duplicateNodes(const int *nodeIdsToDuplicateBg, const int
  *         See \ref MEDCouplingArrayRenumbering for more info on renumbering modes.
  *  \throw If the nodal connectivity of cells is not defined.
  *
+ *  \if ENABLE_EXAMPLES
  *  \ref cpp_mcumesh_renumberNodesInConn "Here is a C++ example".<br>
  *  \ref  py_mcumesh_renumberNodesInConn "Here is a Python example".
+ *  \endif
  */
 void MEDCouplingUMesh::renumberNodesInConn(const int *newNodeNumbersO2N)
 {
@@ -2373,7 +2488,7 @@ void MEDCouplingUMesh::renumberNodesInConn(const int *newNodeNumbersO2N)
  * 
  * \param [in] delta specifies the shift size applied to nodeId in nodal connectivity in \b this.
  */
-void MEDCouplingUMesh::shiftNodeNumbersInConn(int delta) throw(INTERP_KERNEL::Exception)
+void MEDCouplingUMesh::shiftNodeNumbersInConn(int delta)
 {
   checkConnectivityFullyDefined();
   int *conn=getNodalConnectivity()->getPointer();
@@ -2408,7 +2523,7 @@ void MEDCouplingUMesh::shiftNodeNumbersInConn(int delta) throw(INTERP_KERNEL::Ex
  * \param [in] nodeIdsToDuplicateEnd end of node ids (excluded) to be duplicated in connectivity only
  * \param [in] offset the offset applied to all node ids in connectivity that are in [ \a nodeIdsToDuplicateBg, \a nodeIdsToDuplicateEnd ). 
  */
-void MEDCouplingUMesh::duplicateNodesInConn(const int *nodeIdsToDuplicateBg, const int *nodeIdsToDuplicateEnd, int offset) throw(INTERP_KERNEL::Exception)
+void MEDCouplingUMesh::duplicateNodesInConn(const int *nodeIdsToDuplicateBg, const int *nodeIdsToDuplicateEnd, int offset)
 {
   checkConnectivityFullyDefined();
   std::map<int,int> m;
@@ -2448,7 +2563,7 @@ void MEDCouplingUMesh::duplicateNodesInConn(const int *nodeIdsToDuplicateBg, con
  * 
  * \param [in] old2NewBg is expected to be a dynamically allocated pointer of size at least equal to this->getNumberOfCells()
  */
-void MEDCouplingUMesh::renumberCells(const int *old2NewBg, bool check) throw(INTERP_KERNEL::Exception)
+void MEDCouplingUMesh::renumberCells(const int *old2NewBg, bool check)
 {
   checkConnectivityFullyDefined();
   int nbCells=getNumberOfCells();
@@ -2500,8 +2615,10 @@ void MEDCouplingUMesh::renumberCells(const int *old2NewBg, bool check) throw(INT
  *  \throw If the coordinates array is not set.
  *  \throw If the nodal connectivity of cells is not defined.
  *
+ *  \if ENABLE_EXAMPLES
  *  \ref cpp_mcumesh_getCellsInBoundingBox "Here is a C++ example".<br>
  *  \ref  py_mcumesh_getCellsInBoundingBox "Here is a Python example".
+ *  \endif
  */
 DataArrayInt *MEDCouplingUMesh::getCellsInBoundingBox(const double *bbox, double eps) const
 {
@@ -2628,7 +2745,7 @@ INTERP_KERNEL::NormalizedCellType MEDCouplingUMesh::getTypeOfCell(int cellId) co
  * \param [in] type the geometric type
  * \return cell ids in this having geometric type \a type.
  */
-DataArrayInt *MEDCouplingUMesh::giveCellsWithType(INTERP_KERNEL::NormalizedCellType type) const throw(INTERP_KERNEL::Exception)
+DataArrayInt *MEDCouplingUMesh::giveCellsWithType(INTERP_KERNEL::NormalizedCellType type) const
 {
   
   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret=DataArrayInt::New();
@@ -2744,7 +2861,7 @@ std::string MEDCouplingUMesh::advancedRepr() const
  * This method returns a C++ code that is a dump of \a this.
  * This method will throw if this is not fully defined.
  */
-std::string MEDCouplingUMesh::cppRepr() const throw(INTERP_KERNEL::Exception)
+std::string MEDCouplingUMesh::cppRepr() const
 {
   static const char coordsName[]="coords";
   static const char connName[]="conn";
@@ -2778,7 +2895,7 @@ std::string MEDCouplingUMesh::reprConnectivityOfThis() const
  * This method analyzes the 3 arrays of \a this. For each the following behaviour is done : if the array is null a newly one is created
  * with number of tuples set to 0, if not the array is taken as this in the returned instance.
  */
-MEDCouplingUMesh *MEDCouplingUMesh::buildSetInstanceFromThis(int spaceDim) const throw(INTERP_KERNEL::Exception)
+MEDCouplingUMesh *MEDCouplingUMesh::buildSetInstanceFromThis(int spaceDim) const
 {
   int mdim=getMeshDimension();
   if(mdim<0)
@@ -2857,9 +2974,9 @@ int MEDCouplingUMesh::getNumberOfNodesInCell(int cellId) const
  *         describing the cell types. 
  *  \throw If the coordinates array is not set.
  *  \throw If the nodal connectivity of cells is not defined.
- *  \sa getAllTypes()
+ *  \sa getAllGeoTypes()
  */
-std::set<INTERP_KERNEL::NormalizedCellType> MEDCouplingUMesh::getTypesOfPart(const int *begin, const int *end) const throw(INTERP_KERNEL::Exception)
+std::set<INTERP_KERNEL::NormalizedCellType> MEDCouplingUMesh::getTypesOfPart(const int *begin, const int *end) const
 {
   checkFullyDefined();
   std::set<INTERP_KERNEL::NormalizedCellType> ret;
@@ -2932,7 +3049,7 @@ void MEDCouplingUMesh::computeTypes()
 /*!
  * This method checks that all arrays are set. If yes nothing done if no an exception is thrown.
  */
-void MEDCouplingUMesh::checkFullyDefined() const throw(INTERP_KERNEL::Exception)
+void MEDCouplingUMesh::checkFullyDefined() const
 {
   if(!_nodal_connec_index || !_nodal_connec || !_coords)
     throw INTERP_KERNEL::Exception("Reverse nodal connectivity computation requires full connectivity and coordinates set in unstructured mesh.");
@@ -2941,7 +3058,7 @@ void MEDCouplingUMesh::checkFullyDefined() const throw(INTERP_KERNEL::Exception)
 /*!
  * This method checks that all connectivity arrays are set. If yes nothing done if no an exception is thrown.
  */
-void MEDCouplingUMesh::checkConnectivityFullyDefined() const throw(INTERP_KERNEL::Exception)
+void MEDCouplingUMesh::checkConnectivityFullyDefined() const
 {
   if(!_nodal_connec_index || !_nodal_connec)
     throw INTERP_KERNEL::Exception("Reverse nodal connectivity computation requires full connectivity set in unstructured mesh.");
@@ -3113,7 +3230,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;
@@ -3171,7 +3288,7 @@ MEDCouplingFieldDouble *MEDCouplingUMesh::getMeasureField(bool isAbs) const
   name+=getName();
   int nbelem=getNumberOfCells();
   MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> field=MEDCouplingFieldDouble::New(ON_CELLS,ONE_TIME);
-  field->setName(name.c_str());
+  field->setName(name);
   MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> array=DataArrayDouble::New();
   array->alloc(nbelem,1);
   double *area_vol=array->getPointer();
@@ -3216,8 +3333,10 @@ MEDCouplingFieldDouble *MEDCouplingUMesh::getMeasureField(bool isAbs) const
  *  \return DataArrayDouble * - a new instance of DataArrayDouble. The caller is to
  *          delete this array using decrRef() as it is no more needed.
  * 
+ *  \if ENABLE_EXAMPLES
  *  \ref cpp_mcumesh_getPartMeasureField "Here is a C++ example".<br>
  *  \ref  py_mcumesh_getPartMeasureField "Here is a Python example".
+ *  \endif
  *  \sa getMeasureField()
  */
 DataArrayDouble *MEDCouplingUMesh::getPartMeasureField(bool isAbs, const int *begin, const int *end) const
@@ -3226,7 +3345,7 @@ DataArrayDouble *MEDCouplingUMesh::getPartMeasureField(bool isAbs, const int *be
   name+=getName();
   int nbelem=(int)std::distance(begin,end);
   MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> array=DataArrayDouble::New();
-  array->setName(name.c_str());
+  array->setName(name);
   array->alloc(nbelem,1);
   double *area_vol=array->getPointer();
   if(getMeshDimension()!=-1)
@@ -3385,8 +3504,10 @@ MEDCouplingFieldDouble *MEDCouplingUMesh::buildOrthogonalField() const
  *  \throw If the mesh and space dimension is not as specified above.
  *  \sa buildOrthogonalField()
  *
+ *  \if ENABLE_EXAMPLES
  *  \ref cpp_mcumesh_buildPartOrthogonalField "Here is a C++ example".<br>
  *  \ref  py_mcumesh_buildPartOrthogonalField "Here is a Python example".
+ *  \endif
  */
 MEDCouplingFieldDouble *MEDCouplingUMesh::buildPartOrthogonalField(const int *begin, const int *end) const
 {
@@ -3502,7 +3623,7 @@ MEDCouplingFieldDouble *MEDCouplingUMesh::buildDirectionVectorField() const
  *  \throw If the plane does not intersect any 3D cell of \a this mesh.
  *  \throw If \a this includes quadratic cells.
  */
-MEDCouplingUMesh *MEDCouplingUMesh::buildSlice3D(const double *origin, const double *vec, double eps, DataArrayInt *&cellIds) const throw(INTERP_KERNEL::Exception)
+MEDCouplingUMesh *MEDCouplingUMesh::buildSlice3D(const double *origin, const double *vec, double eps, DataArrayInt *&cellIds) const
 {
   checkFullyDefined();
   if(getMeshDimension()!=3 || getSpaceDimension()!=3)
@@ -3568,7 +3689,7 @@ the result mesh.
  *  \throw If the plane does not intersect any 2D cell of \a this mesh.
  *  \throw If \a this includes quadratic cells.
  */
-MEDCouplingUMesh *MEDCouplingUMesh::buildSlice3DSurf(const double *origin, const double *vec, double eps, DataArrayInt *&cellIds) const throw(INTERP_KERNEL::Exception)
+MEDCouplingUMesh *MEDCouplingUMesh::buildSlice3DSurf(const double *origin, const double *vec, double eps, DataArrayInt *&cellIds) const
 {
   checkFullyDefined();
   if(getMeshDimension()!=2 || getSpaceDimension()!=3)
@@ -3649,7 +3770,7 @@ MEDCouplingUMesh *MEDCouplingUMesh::buildSlice3DSurf(const double *origin, const
  *  \throw If magnitude of \a vec is less than 1e-6.
  *  \sa buildSlice3D()
  */
-DataArrayInt *MEDCouplingUMesh::getCellIdsCrossingPlane(const double *origin, const double *vec, double eps) const throw(INTERP_KERNEL::Exception)
+DataArrayInt *MEDCouplingUMesh::getCellIdsCrossingPlane(const double *origin, const double *vec, double eps) const
 {
   checkFullyDefined();
   if(getSpaceDimension()!=3)
@@ -3688,7 +3809,7 @@ DataArrayInt *MEDCouplingUMesh::getCellIdsCrossingPlane(const double *origin, co
  * A 1D mesh is said contiguous if : a cell i with nodal connectivity (k,p) the cell i+1 the nodal connectivity should be (p,m)
  * If not false is returned. In case that false is returned a call to ParaMEDMEM::MEDCouplingUMesh::mergeNodes could be usefull.
  */
-bool MEDCouplingUMesh::isContiguous1D() const throw(INTERP_KERNEL::Exception)
+bool MEDCouplingUMesh::isContiguous1D() const
 {
   if(getMeshDimension()!=1)
     throw INTERP_KERNEL::Exception("MEDCouplingUMesh::isContiguous1D : this method has a sense only for 1D mesh !");
@@ -3767,7 +3888,7 @@ void MEDCouplingUMesh::project1D(const double *pt, const double *v, double eps,
  * dimension - 1.
  * \sa DataArrayDouble::distanceToTuple, MEDCouplingUMesh::distanceToPoints
  */
-double MEDCouplingUMesh::distanceToPoint(const double *ptBg, const double *ptEnd, int& cellId) const throw(INTERP_KERNEL::Exception)
+double MEDCouplingUMesh::distanceToPoint(const double *ptBg, const double *ptEnd, int& cellId) const
 {
   int meshDim=getMeshDimension(),spaceDim=getSpaceDimension();
   if(meshDim!=spaceDim-1)
@@ -3806,7 +3927,7 @@ double MEDCouplingUMesh::distanceToPoint(const double *ptBg, const double *ptEnd
  * \throw if mesh dimension of \a this is not equal to space dimension - 1.
  * \sa DataArrayDouble::distanceToTuple, MEDCouplingUMesh::distanceToPoint
  */
-DataArrayDouble *MEDCouplingUMesh::distanceToPoints(const DataArrayDouble *pts, DataArrayInt *& cellIds) const throw(INTERP_KERNEL::Exception)
+DataArrayDouble *MEDCouplingUMesh::distanceToPoints(const DataArrayDouble *pts, DataArrayInt *& cellIds) const
 {
   if(!pts)
     throw INTERP_KERNEL::Exception("MEDCouplingUMesh::distanceToPoints : input points pointer is NULL !");
@@ -3830,13 +3951,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 +3970,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();
@@ -3877,7 +3998,7 @@ DataArrayDouble *MEDCouplingUMesh::distanceToPoints(const DataArrayDouble *pts,
  * \param [out] cellId that corresponds to minimal distance. If the closer node is not linked to any cell in \a this -1 is returned.
  * \sa MEDCouplingUMesh::distanceToPoint, MEDCouplingUMesh::distanceToPoints
  */
-void MEDCouplingUMesh::DistanceToPoint3DSurfAlg(const double *pt, const int *cellIdsBg, const int *cellIdsEnd, const double *coords, const int *nc, const int *ncI, double& ret0, int& cellId) throw(INTERP_KERNEL::Exception)
+void MEDCouplingUMesh::DistanceToPoint3DSurfAlg(const double *pt, const int *cellIdsBg, const int *cellIdsEnd, const double *coords, const int *nc, const int *ncI, double& ret0, int& cellId)
 {
   cellId=-1;
   ret0=std::numeric_limits<double>::max();
@@ -3916,7 +4037,7 @@ void MEDCouplingUMesh::DistanceToPoint3DSurfAlg(const double *pt, const int *cel
  * \param [out] cellId that corresponds to minimal distance. If the closer node is not linked to any cell in \a this -1 is returned.
  * \sa MEDCouplingUMesh::distanceToPoint, MEDCouplingUMesh::distanceToPoints
  */
-void MEDCouplingUMesh::DistanceToPoint2DCurveAlg(const double *pt, const int *cellIdsBg, const int *cellIdsEnd, const double *coords, const int *nc, const int *ncI, double& ret0, int& cellId) throw(INTERP_KERNEL::Exception)
+void MEDCouplingUMesh::DistanceToPoint2DCurveAlg(const double *pt, const int *cellIdsBg, const int *cellIdsEnd, const double *coords, const int *nc, const int *ncI, double& ret0, int& cellId)
 {
   cellId=-1;
   ret0=std::numeric_limits<double>::max();
@@ -3941,6 +4062,9 @@ void MEDCouplingUMesh::DistanceToPoint2DCurveAlg(const double *pt, const int *ce
 
 /*!
  * Finds cells in contact with a ball (i.e. a point with precision). 
+ * For speed reasons, the INTERP_KERNEL::NORM_QUAD4, INTERP_KERNEL::NORM_TRI6 and INTERP_KERNEL::NORM_QUAD8 cells are considered as convex cells to detect if a point is IN or OUT.
+ * If it is not the case, please change their types to INTERP_KERNEL::NORM_POLYGON or INTERP_KERNEL::NORM_QPOLYG before invoking this method.
+ *
  * \warning This method is suitable if the caller intends to evaluate only one
  *          point, for more points getCellsContainingPoints() is recommended as it is
  *          faster. 
@@ -3962,23 +4086,28 @@ int MEDCouplingUMesh::getCellContainingPoint(const double *pos, double eps) cons
 
 /*!
  * Finds cells in contact with a ball (i.e. a point with precision).
+ * For speed reasons, the INTERP_KERNEL::NORM_QUAD4, INTERP_KERNEL::NORM_TRI6 and INTERP_KERNEL::NORM_QUAD8 cells are considered as convex cells to detect if a point is IN or OUT.
+ * If it is not the case, please change their types to INTERP_KERNEL::NORM_POLYGON or INTERP_KERNEL::NORM_QPOLYG before invoking this method.
  * \warning This method is suitable if the caller intends to evaluate only one
  *          point, for more points getCellsContainingPoints() is recommended as it is
  *          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().
  *
+ *  \if ENABLE_EXAMPLES
  *  \ref cpp_mcumesh_getCellsContainingPoint "Here is a C++ example".<br>
  *  \ref  py_mcumesh_getCellsContainingPoint "Here is a Python example".
+ *  \endif
  */
 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
@@ -4001,6 +4130,41 @@ namespace ParaMEDMEM
     INTERP_KERNEL::NormalizedCellType getTypeOfElement(int) const { return (INTERP_KERNEL::NormalizedCellType)0; }
     // end
   };
+  
+  /*!
+   * Warning the nodes in \a m should be decrRefed ! To avoid that Node * pointer be replaced by another instance.
+   */
+  INTERP_KERNEL::Edge *MEDCouplingUMeshBuildQPFromEdge2(INTERP_KERNEL::NormalizedCellType typ, const int *bg, const double *coords2D, std::map<INTERP_KERNEL::Node *,int>& m)
+  {
+    INTERP_KERNEL::Edge *ret=0;
+    INTERP_KERNEL::Node *n0(new INTERP_KERNEL::Node(coords2D[2*bg[0]],coords2D[2*bg[0]+1])),*n1(new INTERP_KERNEL::Node(coords2D[2*bg[1]],coords2D[2*bg[1]+1]));
+    m[n0]=bg[0]; m[n1]=bg[1];
+    switch(typ)
+      {
+      case INTERP_KERNEL::NORM_SEG2:
+        {
+          ret=new INTERP_KERNEL::EdgeLin(n0,n1);
+          break;
+        }
+      case INTERP_KERNEL::NORM_SEG3:
+        {
+          INTERP_KERNEL::Node *n2(new INTERP_KERNEL::Node(coords2D[2*bg[2]],coords2D[2*bg[2]+1])); m[n2]=bg[2];
+          INTERP_KERNEL::EdgeLin *e1(new INTERP_KERNEL::EdgeLin(n0,n2)),*e2(new INTERP_KERNEL::EdgeLin(n2,n1));
+          INTERP_KERNEL::SegSegIntersector inters(*e1,*e2);
+          // is the SEG3 degenerated, and thus can be reduced to a SEG2?
+          bool colinearity(inters.areColinears());
+          delete e1; delete e2;
+          if(colinearity)
+            { ret=new INTERP_KERNEL::EdgeLin(n0,n1); }
+          else
+            { ret=new INTERP_KERNEL::EdgeArcCircle(n0,n2,n1); }
+          break;
+        }
+      default:
+        throw INTERP_KERNEL::Exception("MEDCouplingUMeshBuildQPFromEdge2 : Expecting a mesh with spaceDim==2 and meshDim==1 !");
+      } 
+    return ret;
+  }
 
   INTERP_KERNEL::Edge *MEDCouplingUMeshBuildQPFromEdge(INTERP_KERNEL::NormalizedCellType typ, std::map<int, std::pair<INTERP_KERNEL::Node *,bool> >& mapp2, const int *bg)
   {
@@ -4017,6 +4181,7 @@ namespace ParaMEDMEM
           INTERP_KERNEL::EdgeLin *e1=new INTERP_KERNEL::EdgeLin(mapp2[bg[0]].first,mapp2[bg[2]].first);
           INTERP_KERNEL::EdgeLin *e2=new INTERP_KERNEL::EdgeLin(mapp2[bg[2]].first,mapp2[bg[1]].first);
           INTERP_KERNEL::SegSegIntersector inters(*e1,*e2);
+          // is the SEG3 degenerated, and thus can be reduced to a SEG2?
           bool colinearity=inters.areColinears();
           delete e1; delete e2;
           if(colinearity)
@@ -4033,11 +4198,14 @@ namespace ParaMEDMEM
   }
 
   /*!
-   * This method creates a sub mesh in Geometric2D DS. The sub mesh is composed be the sub set of cells in 'candidates' and the global mesh 'mDesc'.
-   * The input meth 'mDesc' must be so that mDim==1 et spaceDim==3.
-   * 'mapp' contains a mapping between local numbering in submesh and the global node numbering in 'mDesc'.
+   * This method creates a sub mesh in Geometric2D DS. The sub mesh is composed by the sub set of cells in 'candidates' taken from
+   * the global mesh 'mDesc'.
+   * The input mesh 'mDesc' must be so that mDim==1 and spaceDim==2.
+   * 'mapp' returns a mapping between local numbering in submesh (represented by a Node*) and the global node numbering in 'mDesc'.
    */
-  INTERP_KERNEL::QuadraticPolygon *MEDCouplingUMeshBuildQPFromMesh(const MEDCouplingUMesh *mDesc, const std::vector<int>& candidates, std::map<INTERP_KERNEL::Node *,int>& mapp) throw(INTERP_KERNEL::Exception)
+  INTERP_KERNEL::QuadraticPolygon *MEDCouplingUMeshBuildQPFromMesh(const MEDCouplingUMesh *mDesc, const std::vector<int>& candidates,
+      std::map<INTERP_KERNEL::Node *,int>& mapp)
+      throw(INTERP_KERNEL::Exception)
   {
     mapp.clear();
     std::map<int, std::pair<INTERP_KERNEL::Node *,bool> > mapp2;//bool is for a flag specifying if node is boundary (true) or only a middle for SEG3.
@@ -4082,6 +4250,9 @@ namespace ParaMEDMEM
     return new INTERP_KERNEL::Node(coo1[2*nodeId],coo1[2*nodeId+1]);
   }
 
+  /**
+   * Construct a mapping between set of Nodes and the standart MEDCoupling connectivity format (c, cI).
+   */
   void MEDCouplingUMeshBuildQPFromMesh3(const double *coo1, int offset1, const double *coo2, int offset2, const std::vector<double>& addCoo,
                                         const int *desc1Bg, const int *desc1End, const std::vector<std::vector<int> >& intesctEdges1,
                                         /*output*/std::map<INTERP_KERNEL::Node *,int>& mapp, std::map<int,INTERP_KERNEL::Node *>& mappRev)
@@ -4107,13 +4278,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(eps));
+  const double *bbox(bboxArr->begin());
   int nbOfCells=getNumberOfCells();
   const int *conn=_nodal_connec->getConstPointer();
   const int *connI=_nodal_connec_index->getConstPointer();
@@ -4121,7 +4291,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];
@@ -4131,13 +4301,38 @@ void MEDCouplingUMesh::getCellsContainingPointsAlg(const double *coords, const d
       myTree.getIntersectingElems(bb,candidates);
       for(std::vector<int>::const_iterator iter=candidates.begin();iter!=candidates.end();iter++)
         {
-          int sz=connI[(*iter)+1]-connI[*iter]-1;
-          if(INTERP_KERNEL::PointLocatorAlgos<DummyClsMCUG<SPACEDIM> >::isElementContainsPoint(pos+i*SPACEDIM,
-                                                                                               (INTERP_KERNEL::NormalizedCellType)conn[connI[*iter]],
-                                                                                               coords,conn+connI[*iter]+1,sz,eps))
+          int sz(connI[(*iter)+1]-connI[*iter]-1);
+          INTERP_KERNEL::NormalizedCellType ct((INTERP_KERNEL::NormalizedCellType)conn[connI[*iter]]);
+          bool status(false);
+          if(ct!=INTERP_KERNEL::NORM_POLYGON && ct!=INTERP_KERNEL::NORM_QPOLYG)
+            status=INTERP_KERNEL::PointLocatorAlgos<DummyClsMCUG<SPACEDIM> >::isElementContainsPoint(pos+i*SPACEDIM,ct,coords,conn+connI[*iter]+1,sz,eps);
+          else
+            {
+              if(SPACEDIM!=2)
+                throw INTERP_KERNEL::Exception("MEDCouplingUMesh::getCellsContainingPointsAlg : not implemented yet for POLYGON and QPOLYGON in spaceDim 3 !");
+              INTERP_KERNEL::QUADRATIC_PLANAR::_precision=eps;
+              INTERP_KERNEL::QUADRATIC_PLANAR::_arc_detection_precision=eps;
+              std::vector<INTERP_KERNEL::Node *> nodes(sz);
+              INTERP_KERNEL::QuadraticPolygon *pol(0);
+              for(int j=0;j<sz;j++)
+                {
+                  int nodeId(conn[connI[*iter]+1+j]);
+                  nodes[j]=new INTERP_KERNEL::Node(coords[nodeId*SPACEDIM],coords[nodeId*SPACEDIM+1]);
+                }
+              if(!INTERP_KERNEL::CellModel::GetCellModel(ct).isQuadratic())
+                pol=INTERP_KERNEL::QuadraticPolygon::BuildLinearPolygon(nodes);
+              else
+                pol=INTERP_KERNEL::QuadraticPolygon::BuildArcCirclePolygon(nodes);
+              INTERP_KERNEL::Node *n(new INTERP_KERNEL::Node(pos[i*SPACEDIM],pos[i*SPACEDIM+1]));
+              double a(0.),b(0.),c(0.);
+              a=pol->normalizeMe(b,c); n->applySimilarity(b,c,a);
+              status=pol->isInOrOut2(n);
+              delete pol; n->decrRef();
+            }
+          if(status)
             {
-              eltsIndex[i+1]++;
-              elts.push_back(*iter);
+              eltsIndexPtr[i+1]++;
+              elts->pushBackSilent(*iter);
             }
         }
     }
@@ -4146,13 +4341,15 @@ void MEDCouplingUMesh::getCellsContainingPointsAlg(const double *coords, const d
  * Finds cells in contact with several balls (i.e. points with precision).
  * This method is an extension of getCellContainingPoint() and
  * getCellsContainingPoint() for the case of multiple points.
+ * For speed reasons, the INTERP_KERNEL::NORM_QUAD4, INTERP_KERNEL::NORM_TRI6 and INTERP_KERNEL::NORM_QUAD8 cells are considered as convex cells to detect if a point is IN or OUT.
+ * If it is not the case, please change their types to INTERP_KERNEL::NORM_POLYGON or INTERP_KERNEL::NORM_QPOLYG before invoking this method.
  *  \param [in] pos - an array of coordinates of points in full interlace mode :
  *         X0,Y0,Z0,X1,Y1,Z1,... Size of the array must be \a
  *         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
@@ -4164,11 +4361,13 @@ void MEDCouplingUMesh::getCellsContainingPointsAlg(const double *coords, const d
  *  \throw If the coordinates array is not set.
  *  \throw If \a this->getMeshDimension() != \a this->getSpaceDimension().
  *
+ *  \if ENABLE_EXAMPLES
  *  \ref cpp_mcumesh_getCellsContainingPoints "Here is a C++ example".<br>
  *  \ref  py_mcumesh_getCellsContainingPoints "Here is a Python example".
+ *  \endif
  */
 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();
@@ -4260,7 +4459,7 @@ void MEDCouplingUMesh::checkButterflyCells(std::vector<int>& cells, double eps)
  *
  * \return a newly allocated array containing cellIds that have been modified if any. If no cells have been impacted by this method NULL is returned.
  */
-DataArrayInt *MEDCouplingUMesh::convexEnvelop2D() throw(INTERP_KERNEL::Exception)
+DataArrayInt *MEDCouplingUMesh::convexEnvelop2D()
 {
   if(getMeshDimension()!=2 || getSpaceDimension()!=2)
     throw INTERP_KERNEL::Exception("MEDCouplingUMesh::convexEnvelop2D  works only for meshDim=2 and spaceDim=2 !");
@@ -4355,7 +4554,7 @@ MEDCouplingUMesh *MEDCouplingUMesh::buildExtrudedMesh(const MEDCouplingUMesh *me
  * if cut3DCurve[i]==-1, it means that cell#i has been already detected to be fully part of plane defined by ('origin','vec').
  * This method will throw an exception if \a this contains a non linear segment.
  */
-void MEDCouplingUMesh::split3DCurveWithPlane(const double *origin, const double *vec, double eps, std::vector<int>& cut3DCurve) throw(INTERP_KERNEL::Exception)
+void MEDCouplingUMesh::split3DCurveWithPlane(const double *origin, const double *vec, double eps, std::vector<int>& cut3DCurve)
 {
   checkFullyDefined();
   if(getMeshDimension()!=1 || getSpaceDimension()!=3)
@@ -4462,7 +4661,7 @@ DataArrayDouble *MEDCouplingUMesh::fillExtCoordsUsingTranslation(const MEDCoupli
  * \param mesh1D is the input 1D mesh used for translation and automatic rotation computation.
  * \return newCoords new coords filled by this method. 
  */
-DataArrayDouble *MEDCouplingUMesh::fillExtCoordsUsingTranslAndAutoRotation(const MEDCouplingUMesh *mesh1D, bool isQuad) const throw(INTERP_KERNEL::Exception)
+DataArrayDouble *MEDCouplingUMesh::fillExtCoordsUsingTranslAndAutoRotation(const MEDCouplingUMesh *mesh1D, bool isQuad) const
 {
   if(mesh1D->getSpaceDimension()==2)
     return fillExtCoordsUsingTranslAndAutoRotation2D(mesh1D,isQuad);
@@ -4476,7 +4675,7 @@ DataArrayDouble *MEDCouplingUMesh::fillExtCoordsUsingTranslAndAutoRotation(const
  * \param mesh1D is the input 1D mesh used for translation and automatic rotation computation.
  * \return newCoords new coords filled by this method. 
  */
-DataArrayDouble *MEDCouplingUMesh::fillExtCoordsUsingTranslAndAutoRotation2D(const MEDCouplingUMesh *mesh1D, bool isQuad) const throw(INTERP_KERNEL::Exception)
+DataArrayDouble *MEDCouplingUMesh::fillExtCoordsUsingTranslAndAutoRotation2D(const MEDCouplingUMesh *mesh1D, bool isQuad) const
 {
   if(isQuad)
     throw INTERP_KERNEL::Exception("MEDCouplingUMesh::fillExtCoordsUsingTranslAndAutoRotation2D : not implemented for quadratic cells !");
@@ -4520,7 +4719,7 @@ DataArrayDouble *MEDCouplingUMesh::fillExtCoordsUsingTranslAndAutoRotation2D(con
  * \param mesh1D is the input 1D mesh used for translation and automatic rotation computation.
  * \return newCoords new coords filled by this method. 
  */
-DataArrayDouble *MEDCouplingUMesh::fillExtCoordsUsingTranslAndAutoRotation3D(const MEDCouplingUMesh *mesh1D, bool isQuad) const throw(INTERP_KERNEL::Exception)
+DataArrayDouble *MEDCouplingUMesh::fillExtCoordsUsingTranslAndAutoRotation3D(const MEDCouplingUMesh *mesh1D, bool isQuad) const
 {
   if(isQuad)
     throw INTERP_KERNEL::Exception("MEDCouplingUMesh::fillExtCoordsUsingTranslAndAutoRotation3D : not implemented for quadratic cells !");
@@ -4680,7 +4879,7 @@ bool MEDCouplingUMesh::isPresenceOfQuadratic() const
  *  \throw If the coordinates array is not set.
  *  \throw If the nodal connectivity of cells is not defined.
  */
-void MEDCouplingUMesh::convertQuadraticCellsToLinear() throw(INTERP_KERNEL::Exception)
+void MEDCouplingUMesh::convertQuadraticCellsToLinear()
 {
   checkFullyDefined();
   int nbOfCells=getNumberOfCells();
@@ -4753,7 +4952,7 @@ void MEDCouplingUMesh::convertQuadraticCellsToLinear() throw(INTERP_KERNEL::Exce
  *
  * \sa MEDCouplingUMesh::convertQuadraticCellsToLinear
  */
-DataArrayInt *MEDCouplingUMesh::convertLinearCellsToQuadratic(int conversionType) throw(INTERP_KERNEL::Exception)
+DataArrayInt *MEDCouplingUMesh::convertLinearCellsToQuadratic(int conversionType)
 {
   DataArrayInt *conn=0,*connI=0;
   DataArrayDouble *coords=0;
@@ -4818,7 +5017,7 @@ DataArrayInt *MEDCouplingUMesh::convertLinearCellsToQuadratic(int conversionType
  * \return a newly created DataArrayInt instance that the caller should deal with containing cell ids of converted cells.
  * \sa MEDCouplingUMesh::convertLinearCellsToQuadratic.
  */
-DataArrayInt *MEDCouplingUMesh::convertLinearCellsToQuadratic1D0(DataArrayInt *&conn, DataArrayInt *&connI, DataArrayDouble *& coords, std::set<INTERP_KERNEL::NormalizedCellType>& types) const throw(INTERP_KERNEL::Exception)
+DataArrayInt *MEDCouplingUMesh::convertLinearCellsToQuadratic1D0(DataArrayInt *&conn, DataArrayInt *&connI, DataArrayDouble *& coords, std::set<INTERP_KERNEL::NormalizedCellType>& types) const
 {
   MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> bary=getBarycenterAndOwner();
   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> newConn=DataArrayInt::New(); newConn->alloc(0,1);
@@ -4855,7 +5054,7 @@ DataArrayInt *MEDCouplingUMesh::convertLinearCellsToQuadratic1D0(DataArrayInt *&
   return ret.retn();
 }
 
-DataArrayInt *MEDCouplingUMesh::convertLinearCellsToQuadratic2DAnd3D0(const MEDCouplingUMesh *m1D, const DataArrayInt *desc, const DataArrayInt *descI, DataArrayInt *&conn, DataArrayInt *&connI, DataArrayDouble *& coords, std::set<INTERP_KERNEL::NormalizedCellType>& types) const throw(INTERP_KERNEL::Exception)
+DataArrayInt *MEDCouplingUMesh::convertLinearCellsToQuadratic2DAnd3D0(const MEDCouplingUMesh *m1D, const DataArrayInt *desc, const DataArrayInt *descI, DataArrayInt *&conn, DataArrayInt *&connI, DataArrayDouble *& coords, std::set<INTERP_KERNEL::NormalizedCellType>& types) const
 {
   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> newConn=DataArrayInt::New(); newConn->alloc(0,1);
   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> newConnI=DataArrayInt::New(); newConnI->alloc(1,1); newConnI->setIJ(0,0,0);
@@ -4906,7 +5105,7 @@ DataArrayInt *MEDCouplingUMesh::convertLinearCellsToQuadratic2DAnd3D0(const MEDC
  * \return a newly created DataArrayInt instance that the caller should deal with containing cell ids of converted cells.
  * \sa MEDCouplingUMesh::convertLinearCellsToQuadratic.
  */
-DataArrayInt *MEDCouplingUMesh::convertLinearCellsToQuadratic2D0(DataArrayInt *&conn, DataArrayInt *&connI, DataArrayDouble *& coords, std::set<INTERP_KERNEL::NormalizedCellType>& types) const throw(INTERP_KERNEL::Exception)
+DataArrayInt *MEDCouplingUMesh::convertLinearCellsToQuadratic2D0(DataArrayInt *&conn, DataArrayInt *&connI, DataArrayDouble *& coords, std::set<INTERP_KERNEL::NormalizedCellType>& types) const
 {
   
   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> desc(DataArrayInt::New()),descI(DataArrayInt::New()),tmp2(DataArrayInt::New()),tmp3(DataArrayInt::New());
@@ -4914,7 +5113,7 @@ DataArrayInt *MEDCouplingUMesh::convertLinearCellsToQuadratic2D0(DataArrayInt *&
   return convertLinearCellsToQuadratic2DAnd3D0(m1D,desc,descI,conn,connI,coords,types);
 }
 
-DataArrayInt *MEDCouplingUMesh::convertLinearCellsToQuadratic2D1(DataArrayInt *&conn, DataArrayInt *&connI, DataArrayDouble *& coords, std::set<INTERP_KERNEL::NormalizedCellType>& types) const throw(INTERP_KERNEL::Exception)
+DataArrayInt *MEDCouplingUMesh::convertLinearCellsToQuadratic2D1(DataArrayInt *&conn, DataArrayInt *&connI, DataArrayDouble *& coords, std::set<INTERP_KERNEL::NormalizedCellType>& types) const
 {
   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> desc(DataArrayInt::New()),descI(DataArrayInt::New()),tmp2(DataArrayInt::New()),tmp3(DataArrayInt::New());
   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> m1D=buildDescendingConnectivity(desc,descI,tmp2,tmp3); tmp2=0; tmp3=0;
@@ -4971,14 +5170,14 @@ DataArrayInt *MEDCouplingUMesh::convertLinearCellsToQuadratic2D1(DataArrayInt *&
  * \return a newly created DataArrayInt instance that the caller should deal with containing cell ids of converted cells.
  * \sa MEDCouplingUMesh::convertLinearCellsToQuadratic.
  */
-DataArrayInt *MEDCouplingUMesh::convertLinearCellsToQuadratic3D0(DataArrayInt *&conn, DataArrayInt *&connI, DataArrayDouble *& coords, std::set<INTERP_KERNEL::NormalizedCellType>& types) const throw(INTERP_KERNEL::Exception)
+DataArrayInt *MEDCouplingUMesh::convertLinearCellsToQuadratic3D0(DataArrayInt *&conn, DataArrayInt *&connI, DataArrayDouble *& coords, std::set<INTERP_KERNEL::NormalizedCellType>& types) const
 {
   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> desc(DataArrayInt::New()),descI(DataArrayInt::New()),tmp2(DataArrayInt::New()),tmp3(DataArrayInt::New());
   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> m1D=explode3DMeshTo1D(desc,descI,tmp2,tmp3); tmp2=0; tmp3=0;
   return convertLinearCellsToQuadratic2DAnd3D0(m1D,desc,descI,conn,connI,coords,types);
 }
 
-DataArrayInt *MEDCouplingUMesh::convertLinearCellsToQuadratic3D1(DataArrayInt *&conn, DataArrayInt *&connI, DataArrayDouble *& coords, std::set<INTERP_KERNEL::NormalizedCellType>& types) const throw(INTERP_KERNEL::Exception)
+DataArrayInt *MEDCouplingUMesh::convertLinearCellsToQuadratic3D1(DataArrayInt *&conn, DataArrayInt *&connI, DataArrayDouble *& coords, std::set<INTERP_KERNEL::NormalizedCellType>& types) const
 {
   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> desc2(DataArrayInt::New()),desc2I(DataArrayInt::New()),tmp2(DataArrayInt::New()),tmp3(DataArrayInt::New());
   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> m2D=buildDescendingConnectivityGen<MinusOneSonsGeneratorBiQuadratic>(desc2,desc2I,tmp2,tmp3,MEDCouplingFastNbrer); tmp2=0; tmp3=0;
@@ -5071,7 +5270,7 @@ DataArrayInt *MEDCouplingUMesh::convertLinearCellsToQuadratic3D1(DataArrayInt *&
  *  \throw If \a this->getMeshDimension() != 2.
  *  \throw If \a this->getSpaceDimension() != 2.
  */
-void MEDCouplingUMesh::tessellate2D(double eps) throw(INTERP_KERNEL::Exception)
+void MEDCouplingUMesh::tessellate2D(double eps)
 {
   checkFullyDefined();
   if(getMeshDimension()!=2 || getSpaceDimension()!=2)  
@@ -5100,7 +5299,7 @@ void MEDCouplingUMesh::tessellate2D(double eps) throw(INTERP_KERNEL::Exception)
  *  \throw If \a this->getMeshDimension() != 1.
  *  \throw If \a this->getSpaceDimension() != 2.
  */
-void MEDCouplingUMesh::tessellate2DCurve(double eps) throw(INTERP_KERNEL::Exception)
+void MEDCouplingUMesh::tessellate2DCurve(double eps)
 {
   checkFullyDefined();
   if(getMeshDimension()!=1 || getSpaceDimension()!=2)
@@ -5155,7 +5354,7 @@ void MEDCouplingUMesh::tessellate2DCurve(double eps) throw(INTERP_KERNEL::Except
           newConnIPtr[1]=newConnIPtr[0]+3;
         }
     }
-  if(addCoo.empty() && ((int)newConn.size())==_nodal_connec->getNumberOfTuples())//nothing happens during tasselation : no update needed
+  if(addCoo.empty() && ((int)newConn.size())==_nodal_connec->getNumberOfTuples())//nothing happens during tessellation : no update needed
     return ;
   _types=types;
   DataArrayInt::SetArrayIn(newConnI,_nodal_connec_index);
@@ -5192,8 +5391,9 @@ 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, MEDCoupling1SGTUMesh::sortHexa8EachOther
  */
-DataArrayInt *MEDCouplingUMesh::simplexize(int policy) throw(INTERP_KERNEL::Exception)
+DataArrayInt *MEDCouplingUMesh::simplexize(int policy)
 {
   switch(policy)
     {
@@ -5224,7 +5424,7 @@ DataArrayInt *MEDCouplingUMesh::simplexize(int policy) throw(INTERP_KERNEL::Exce
  *  \throw If the nodal connectivity of cells is not defined.
  *  \throw If \a this->getMeshDimension() < 1.
  */
-bool MEDCouplingUMesh::areOnlySimplexCells() const throw(INTERP_KERNEL::Exception)
+bool MEDCouplingUMesh::areOnlySimplexCells() const
 {
   checkFullyDefined();
   int mdim=getMeshDimension();
@@ -5245,7 +5445,7 @@ bool MEDCouplingUMesh::areOnlySimplexCells() const throw(INTERP_KERNEL::Exceptio
 /*!
  * This method implements policy 0 of virtual method ParaMEDMEM::MEDCouplingUMesh::simplexize.
  */
-DataArrayInt *MEDCouplingUMesh::simplexizePol0() throw(INTERP_KERNEL::Exception)
+DataArrayInt *MEDCouplingUMesh::simplexizePol0()
 {
   checkConnectivityFullyDefined();
   if(getMeshDimension()!=2)
@@ -5298,7 +5498,7 @@ DataArrayInt *MEDCouplingUMesh::simplexizePol0() throw(INTERP_KERNEL::Exception)
 /*!
  * This method implements policy 1 of virtual method ParaMEDMEM::MEDCouplingUMesh::simplexize.
  */
-DataArrayInt *MEDCouplingUMesh::simplexizePol1() throw(INTERP_KERNEL::Exception)
+DataArrayInt *MEDCouplingUMesh::simplexizePol1()
 {
   checkConnectivityFullyDefined();
   if(getMeshDimension()!=2)
@@ -5351,7 +5551,7 @@ DataArrayInt *MEDCouplingUMesh::simplexizePol1() throw(INTERP_KERNEL::Exception)
 /*!
  * This method implements policy INTERP_KERNEL::PLANAR_FACE_5 of virtual method ParaMEDMEM::MEDCouplingUMesh::simplexize.
  */
-DataArrayInt *MEDCouplingUMesh::simplexizePlanarFace5() throw(INTERP_KERNEL::Exception)
+DataArrayInt *MEDCouplingUMesh::simplexizePlanarFace5()
 {
   checkConnectivityFullyDefined();
   if(getMeshDimension()!=3)
@@ -5403,7 +5603,7 @@ DataArrayInt *MEDCouplingUMesh::simplexizePlanarFace5() throw(INTERP_KERNEL::Exc
 /*!
  * This method implements policy INTERP_KERNEL::PLANAR_FACE_6 of virtual method ParaMEDMEM::MEDCouplingUMesh::simplexize.
  */
-DataArrayInt *MEDCouplingUMesh::simplexizePlanarFace6() throw(INTERP_KERNEL::Exception)
+DataArrayInt *MEDCouplingUMesh::simplexizePlanarFace6()
 {
   checkConnectivityFullyDefined();
   if(getMeshDimension()!=3)
@@ -5460,7 +5660,7 @@ DataArrayInt *MEDCouplingUMesh::simplexizePlanarFace6() throw(INTERP_KERNEL::Exc
  * \param desc is descending connectivity in format specified in MEDCouplingUMesh::buildDescendingConnectivity2
  * \param descIndex is descending connectivity index in format specified in MEDCouplingUMesh::buildDescendingConnectivity2
  */
-void MEDCouplingUMesh::subDivide2DMesh(const int *nodeSubdived, const int *nodeIndxSubdived, const int *desc, const int *descIndex) throw(INTERP_KERNEL::Exception)
+void MEDCouplingUMesh::subDivide2DMesh(const int *nodeSubdived, const int *nodeIndxSubdived, const int *desc, const int *descIndex)
 {
   checkFullyDefined();
   if(getMeshDimension()!=2)
@@ -5544,7 +5744,7 @@ void MEDCouplingUMesh::subDivide2DMesh(const int *nodeSubdived, const int *nodeI
  *  \throw If the coordinates array is not set.
  *  \throw If the nodal connectivity of cells is not defined.
  */
-void MEDCouplingUMesh::convertDegeneratedCells() throw(INTERP_KERNEL::Exception)
+void MEDCouplingUMesh::convertDegeneratedCells()
 {
   checkFullyDefined();
   if(getMeshDimension()<=1)
@@ -5588,10 +5788,12 @@ void MEDCouplingUMesh::convertDegeneratedCells() throw(INTERP_KERNEL::Exception)
  *  \throw If \a this->getMeshDimension() != 2.
  *  \throw If \a this->getSpaceDimension() != 3.
  *
+ *  \if ENABLE_EXAMPLES
  *  \ref cpp_mcumesh_are2DCellsNotCorrectlyOriented "Here is a C++ example".<br>
  *  \ref  py_mcumesh_are2DCellsNotCorrectlyOriented "Here is a Python example".
+ *  \endif
  */
-void MEDCouplingUMesh::are2DCellsNotCorrectlyOriented(const double *vec, bool polyOnly, std::vector<int>& cells) const throw(INTERP_KERNEL::Exception)
+void MEDCouplingUMesh::are2DCellsNotCorrectlyOriented(const double *vec, bool polyOnly, std::vector<int>& cells) const
 {
   if(getMeshDimension()!=2 || getSpaceDimension()!=3)
     throw INTERP_KERNEL::Exception("Invalid mesh to apply are2DCellsNotCorrectlyOriented on it : must be meshDim==2 and spaceDim==3 !");
@@ -5622,10 +5824,12 @@ void MEDCouplingUMesh::are2DCellsNotCorrectlyOriented(const double *vec, bool po
  *  \throw If \a this->getMeshDimension() != 2.
  *  \throw If \a this->getSpaceDimension() != 3.
  *
+ *  \if ENABLE_EXAMPLES
  *  \ref cpp_mcumesh_are2DCellsNotCorrectlyOriented "Here is a C++ example".<br>
  *  \ref  py_mcumesh_are2DCellsNotCorrectlyOriented "Here is a Python example".
+ *  \endif
  */
-void MEDCouplingUMesh::orientCorrectly2DCells(const double *vec, bool polyOnly) throw(INTERP_KERNEL::Exception)
+void MEDCouplingUMesh::orientCorrectly2DCells(const double *vec, bool polyOnly)
 {
   if(getMeshDimension()!=2 || getSpaceDimension()!=3)
     throw INTERP_KERNEL::Exception("Invalid mesh to apply orientCorrectly2DCells on it : must be meshDim==2 and spaceDim==3 !");
@@ -5639,13 +5843,25 @@ void MEDCouplingUMesh::orientCorrectly2DCells(const double *vec, bool polyOnly)
       INTERP_KERNEL::NormalizedCellType type=(INTERP_KERNEL::NormalizedCellType)conn[connI[i]];
       if(!polyOnly || (type==INTERP_KERNEL::NORM_POLYGON || type==INTERP_KERNEL::NORM_QPOLYG))
         {
-          bool isQuadratic=INTERP_KERNEL::CellModel::GetCellModel(type).isQuadratic();
+          bool isQuadratic(INTERP_KERNEL::CellModel::GetCellModel(type).isQuadratic());
           if(!IsPolygonWellOriented(isQuadratic,vec,conn+connI[i]+1,conn+connI[i+1],coordsPtr))
             {
               isModified=true;
-              std::vector<int> tmp(connI[i+1]-connI[i]-2);
-              std::copy(conn+connI[i]+2,conn+connI[i+1],tmp.rbegin());
-              std::copy(tmp.begin(),tmp.end(),conn+connI[i]+2);
+              if(!isQuadratic)
+                {
+                  std::vector<int> tmp(connI[i+1]-connI[i]-2);
+                  std::copy(conn+connI[i]+2,conn+connI[i+1],tmp.rbegin());
+                  std::copy(tmp.begin(),tmp.end(),conn+connI[i]+2);
+                }
+              else
+                {
+                  int sz(((int)(connI[i+1]-connI[i]-1))/2);
+                  std::vector<int> tmp0(sz-1),tmp1(sz);
+                  std::copy(conn+connI[i]+2,conn+connI[i]+1+sz,tmp0.rbegin());
+                  std::copy(conn+connI[i]+1+sz,conn+connI[i+1],tmp1.rbegin());
+                  std::copy(tmp0.begin(),tmp0.end(),conn+connI[i]+2);
+                  std::copy(tmp1.begin(),tmp1.end(),conn+connI[i]+1+sz);
+                }
             }
         }
     }
@@ -5664,10 +5880,12 @@ void MEDCouplingUMesh::orientCorrectly2DCells(const double *vec, bool polyOnly)
  *  \throw If the coordinates array is not set.
  *  \throw If the nodal connectivity of cells is not defined.
  *
+ *  \if ENABLE_EXAMPLES
  *  \ref cpp_mcumesh_arePolyhedronsNotCorrectlyOriented "Here is a C++ example".<br>
  *  \ref  py_mcumesh_arePolyhedronsNotCorrectlyOriented "Here is a Python example".
+ *  \endif
  */
-void MEDCouplingUMesh::arePolyhedronsNotCorrectlyOriented(std::vector<int>& cells) const throw(INTERP_KERNEL::Exception)
+void MEDCouplingUMesh::arePolyhedronsNotCorrectlyOriented(std::vector<int>& cells) const
 {
   if(getMeshDimension()!=3 || getSpaceDimension()!=3)
     throw INTERP_KERNEL::Exception("Invalid mesh to apply arePolyhedronsNotCorrectlyOriented on it : must be meshDim==3 and spaceDim==3 !");
@@ -5695,11 +5913,13 @@ void MEDCouplingUMesh::arePolyhedronsNotCorrectlyOriented(std::vector<int>& cell
  *  \throw If the nodal connectivity of cells is not defined.
  *  \throw If the reparation fails.
  *
+ *  \if ENABLE_EXAMPLES
  *  \ref cpp_mcumesh_arePolyhedronsNotCorrectlyOriented "Here is a C++ example".<br>
  *  \ref  py_mcumesh_arePolyhedronsNotCorrectlyOriented "Here is a Python example".
+ *  \endif
  * \sa MEDCouplingUMesh::findAndCorrectBadOriented3DCells
  */
-void MEDCouplingUMesh::orientCorrectlyPolyhedrons() throw(INTERP_KERNEL::Exception)
+void MEDCouplingUMesh::orientCorrectlyPolyhedrons()
 {
   if(getMeshDimension()!=3 || getSpaceDimension()!=3)
     throw INTERP_KERNEL::Exception("Invalid mesh to apply orientCorrectlyPolyhedrons on it : must be meshDim==3 and spaceDim==3 !");
@@ -5740,11 +5960,13 @@ void MEDCouplingUMesh::orientCorrectlyPolyhedrons() throw(INTERP_KERNEL::Excepti
  *  \throw If the coordinates array is not set.
  *  \throw If the nodal connectivity of cells is not defined.
  *
+ *  \if ENABLE_EXAMPLES
  *  \ref cpp_mcumesh_findAndCorrectBadOriented3DExtrudedCells "Here is a C++ example".<br>
  *  \ref  py_mcumesh_findAndCorrectBadOriented3DExtrudedCells "Here is a Python example".
+ *  \endif
  * \sa MEDCouplingUMesh::findAndCorrectBadOriented3DCells
  */
-DataArrayInt *MEDCouplingUMesh::findAndCorrectBadOriented3DExtrudedCells() throw(INTERP_KERNEL::Exception)
+DataArrayInt *MEDCouplingUMesh::findAndCorrectBadOriented3DExtrudedCells()
 {
   const char msg[]="check3DCellsWellOriented detection works only for 3D cells !";
   if(getMeshDimension()!=3)
@@ -5781,7 +6003,7 @@ DataArrayInt *MEDCouplingUMesh::findAndCorrectBadOriented3DExtrudedCells() throw
  * \ret a newly allocated int array with one components containing cell ids renumbered to fit the convention of MED (MED file and MEDCoupling)
  * \sa MEDCouplingUMesh::orientCorrectlyPolyhedrons, 
  */
-DataArrayInt *MEDCouplingUMesh::findAndCorrectBadOriented3DCells() throw(INTERP_KERNEL::Exception)
+DataArrayInt *MEDCouplingUMesh::findAndCorrectBadOriented3DCells()
 {
   if(getMeshDimension()!=3 || getSpaceDimension()!=3)
     throw INTERP_KERNEL::Exception("Invalid mesh to apply findAndCorrectBadOriented3DCells on it : must be meshDim==3 and spaceDim==3 !");
@@ -5848,7 +6070,7 @@ DataArrayInt *MEDCouplingUMesh::findAndCorrectBadOriented3DCells() throw(INTERP_
  * \param vec output of size at least 3 used to store the normal vector (with norm equal to Area ) of searched plane.
  * \param pos output of size at least 3 used to store a point owned of searched plane.
  */
-void MEDCouplingUMesh::getFastAveragePlaneOfThis(double *vec, double *pos) const throw(INTERP_KERNEL::Exception)
+void MEDCouplingUMesh::getFastAveragePlaneOfThis(double *vec, double *pos) const
 {
   if(getMeshDimension()!=2 || getSpaceDimension()!=3)
     throw INTERP_KERNEL::Exception("Invalid mesh to apply getFastAveragePlaneOfThis on it : must be meshDim==2 and spaceDim==3 !");
@@ -5882,7 +6104,7 @@ void MEDCouplingUMesh::getFastAveragePlaneOfThis(double *vec, double *pos) const
  *  \throw If \a this->getSpaceDimension() is neither 2 nor 3.
  *  \throw If \a this mesh includes cells of type different from the ones enumerated above.
  */
-MEDCouplingFieldDouble *MEDCouplingUMesh::getEdgeRatioField() const throw(INTERP_KERNEL::Exception)
+MEDCouplingFieldDouble *MEDCouplingUMesh::getEdgeRatioField() const
 {
   checkCoherency();
   int spaceDim=getSpaceDimension();
@@ -5954,7 +6176,7 @@ MEDCouplingFieldDouble *MEDCouplingUMesh::getEdgeRatioField() const throw(INTERP
  *  \throw If \a this->getSpaceDimension() is neither 2 nor 3.
  *  \throw If \a this mesh includes cells of type different from the ones enumerated above.
  */
-MEDCouplingFieldDouble *MEDCouplingUMesh::getAspectRatioField() const throw(INTERP_KERNEL::Exception)
+MEDCouplingFieldDouble *MEDCouplingUMesh::getAspectRatioField() const
 {
   checkCoherency();
   int spaceDim=getSpaceDimension();
@@ -6025,7 +6247,7 @@ MEDCouplingFieldDouble *MEDCouplingUMesh::getAspectRatioField() const throw(INTE
  *  \throw If \a this->getSpaceDimension() != 3.
  *  \throw If \a this mesh includes cells of type different from the ones enumerated above.
  */
-MEDCouplingFieldDouble *MEDCouplingUMesh::getWarpField() const throw(INTERP_KERNEL::Exception)
+MEDCouplingFieldDouble *MEDCouplingUMesh::getWarpField() const
 {
   checkCoherency();
   int spaceDim=getSpaceDimension();
@@ -6085,7 +6307,7 @@ MEDCouplingFieldDouble *MEDCouplingUMesh::getWarpField() const throw(INTERP_KERN
  *  \throw If \a this->getSpaceDimension() != 3.
  *  \throw If \a this mesh includes cells of type different from the ones enumerated above.
  */
-MEDCouplingFieldDouble *MEDCouplingUMesh::getSkewField() const throw(INTERP_KERNEL::Exception)
+MEDCouplingFieldDouble *MEDCouplingUMesh::getSkewField() const
 {
   checkCoherency();
   int spaceDim=getSpaceDimension();
@@ -6128,36 +6350,175 @@ 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.
+ * 
+ * \param [in] arcDetEps - a parameter specifying in case of 2D quadratic polygon cell the detection limit between linear and arc circle. (By default 1e-12)
+ *                         For all other cases this input parameter is ignored.
+ * \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.
+ * \sa MEDCouplingUMesh::getBoundingBoxForBBTreeFast, MEDCouplingUMesh::getBoundingBoxForBBTree2DQuadratic
  */
-void MEDCouplingUMesh::getBoundingBoxForBBTree(std::vector<double>& bbox) const
+DataArrayDouble *MEDCouplingUMesh::getBoundingBoxForBBTree(double arcDetEps) const
 {
-  int spaceDim=getSpaceDimension();
-  int nbOfCells=getNumberOfCells();
-  bbox.resize(2*nbOfCells*spaceDim);
+  int mDim(getMeshDimension()),sDim(getSpaceDimension());
+  if((mDim==3 && sDim==3) || (mDim==2 && sDim==3) || (mDim==1 && sDim==1) || ( mDim==1 && sDim==3))  // Compute refined boundary box for quadratic elements only in 2D.
+    return getBoundingBoxForBBTreeFast();
+  if((mDim==2 && sDim==2) || (mDim==1 && sDim==2))
+    {
+      bool presenceOfQuadratic(false);
+      for(std::set<INTERP_KERNEL::NormalizedCellType>::const_iterator it=_types.begin();it!=_types.end();it++)
+        {
+          const INTERP_KERNEL::CellModel& cm(INTERP_KERNEL::CellModel::GetCellModel(*it));
+          if(cm.isQuadratic())
+            presenceOfQuadratic=true;
+        }
+      if(!presenceOfQuadratic)
+        return getBoundingBoxForBBTreeFast();
+      if(mDim==2 && sDim==2)
+        return getBoundingBoxForBBTree2DQuadratic(arcDetEps);
+      else
+        return getBoundingBoxForBBTree1DQuadratic(arcDetEps);
+    }
+  throw INTERP_KERNEL::Exception("MEDCouplingUMesh::getBoundingBoxForBBTree : Managed dimensions are (mDim=1,sDim=1), (mDim=1,sDim=2), (mDim=1,sDim=3), (mDim=2,sDim=2), (mDim=2,sDim=3) and (mDim=3,sDim=3) !");
+}
+
+/*!
+ * This method aggregate the bbox of each cell only considering the nodes constituting each cell and put it into bbox parameter.
+ * So meshes having quadratic cells the computed bounding boxes can be invalid !
+ * 
+ * \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.
+ */
+DataArrayDouble *MEDCouplingUMesh::getBoundingBoxForBBTreeFast() const
+{
+  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();
+}
+
+/*!
+ * This method aggregates the bbox of each 2D cell in \a this considering the whole shape. This method is particularly
+ * useful for 2D meshes having quadratic cells
+ * because for this type of cells getBoundingBoxForBBTreeFast method may return invalid bounding boxes (since it just considers
+ * the two extremities of the arc of circle).
+ * 
+ * \param [in] arcDetEps - a parameter specifying in case of 2D quadratic polygon cell the detection limit between linear and arc circle. (By default 1e-12)
+ * \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 defined.
+ * \throw If \a this is not a mesh with meshDimension equal to 2.
+ * \throw If \a this is not a mesh with spaceDimension equal to 2.
+ * \sa MEDCouplingUMesh::getBoundingBoxForBBTree1DQuadratic
+ */
+DataArrayDouble *MEDCouplingUMesh::getBoundingBoxForBBTree2DQuadratic(double arcDetEps) const
+{
+  checkFullyDefined();
+  int spaceDim(getSpaceDimension()),mDim(getMeshDimension()),nbOfCells(getNumberOfCells());
+  if(spaceDim!=2 || mDim!=2)
+    throw INTERP_KERNEL::Exception("MEDCouplingUMesh::getBoundingBoxForBBTree2DQuadratic : This method should be applied on mesh with mesh dimension equal to 2 and space dimension also equal to 2!");
+  MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> ret(DataArrayDouble::New()); ret->alloc(nbOfCells,2*spaceDim);
+  double *bbox(ret->getPointer());
+  const double *coords(_coords->getConstPointer());
+  const int *conn(_nodal_connec->getConstPointer()),*connI(_nodal_connec_index->getConstPointer());
+  for(int i=0;i<nbOfCells;i++,bbox+=4,connI++)
+    {
+      const INTERP_KERNEL::CellModel& cm(INTERP_KERNEL::CellModel::GetCellModel((INTERP_KERNEL::NormalizedCellType)conn[*connI]));
+      int sz(connI[1]-connI[0]-1);
+      INTERP_KERNEL::QUADRATIC_PLANAR::_arc_detection_precision=arcDetEps;
+      std::vector<INTERP_KERNEL::Node *> nodes(sz);
+      INTERP_KERNEL::QuadraticPolygon *pol(0);
+      for(int j=0;j<sz;j++)
+        {
+          int nodeId(conn[*connI+1+j]);
+          nodes[j]=new INTERP_KERNEL::Node(coords[nodeId*2],coords[nodeId*2+1]);
+        }
+      if(!cm.isQuadratic())
+        pol=INTERP_KERNEL::QuadraticPolygon::BuildLinearPolygon(nodes);
+      else
+        pol=INTERP_KERNEL::QuadraticPolygon::BuildArcCirclePolygon(nodes);
+      INTERP_KERNEL::Bounds b; pol->fillBounds(b); delete pol;
+      bbox[0]=b.getXMin(); bbox[1]=b.getXMax(); bbox[2]=b.getYMin(); bbox[3]=b.getYMax(); 
     }
+  return ret.retn();
+}
+
+/*!
+ * This method aggregates the bbox of each 1D cell in \a this considering the whole shape. This method is particularly
+ * useful for 2D meshes having quadratic cells
+ * because for this type of cells getBoundingBoxForBBTreeFast method may return invalid bounding boxes (since it just considers
+ * the two extremities of the arc of circle).
+ * 
+ * \param [in] arcDetEps - a parameter specifying in case of 2D quadratic polygon cell the detection limit between linear and arc circle. (By default 1e-12)
+ * \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 defined.
+ * \throw If \a this is not a mesh with meshDimension equal to 1.
+ * \throw If \a this is not a mesh with spaceDimension equal to 2.
+ * \sa MEDCouplingUMesh::getBoundingBoxForBBTree2DQuadratic
+ */
+DataArrayDouble *MEDCouplingUMesh::getBoundingBoxForBBTree1DQuadratic(double arcDetEps) const
+{
+  checkFullyDefined();
+  int spaceDim(getSpaceDimension()),mDim(getMeshDimension()),nbOfCells(getNumberOfCells());
+  if(spaceDim!=2 || mDim!=1)
+    throw INTERP_KERNEL::Exception("MEDCouplingUMesh::getBoundingBoxForBBTree1DQuadratic : This method should be applied on mesh with mesh dimension equal to 1 and space dimension also equal to 2!");
+  MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> ret(DataArrayDouble::New()); ret->alloc(nbOfCells,2*spaceDim);
+  double *bbox(ret->getPointer());
+  const double *coords(_coords->getConstPointer());
+  const int *conn(_nodal_connec->getConstPointer()),*connI(_nodal_connec_index->getConstPointer());
+  for(int i=0;i<nbOfCells;i++,bbox+=4,connI++)
+    {
+      const INTERP_KERNEL::CellModel& cm(INTERP_KERNEL::CellModel::GetCellModel((INTERP_KERNEL::NormalizedCellType)conn[*connI]));
+      int sz(connI[1]-connI[0]-1);
+      INTERP_KERNEL::QUADRATIC_PLANAR::_arc_detection_precision=arcDetEps;
+      std::vector<INTERP_KERNEL::Node *> nodes(sz);
+      INTERP_KERNEL::Edge *edge(0);
+      for(int j=0;j<sz;j++)
+        {
+          int nodeId(conn[*connI+1+j]);
+          nodes[j]=new INTERP_KERNEL::Node(coords[nodeId*2],coords[nodeId*2+1]);
+        }
+      if(!cm.isQuadratic())
+        edge=INTERP_KERNEL::QuadraticPolygon::BuildLinearEdge(nodes);
+      else
+        edge=INTERP_KERNEL::QuadraticPolygon::BuildArcCircleEdge(nodes);
+      const INTERP_KERNEL::Bounds& b(edge->getBounds());
+      bbox[0]=b.getXMin(); bbox[1]=b.getXMax(); bbox[2]=b.getYMin(); bbox[3]=b.getYMax(); edge->decrRef();
+    }
+  return ret.retn();
 }
 
 /// @cond INTERNAL
@@ -6192,18 +6553,18 @@ 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)
+std::vector<int> MEDCouplingUMesh::getDistributionOfTypes() const
 {
   checkConnectivityFullyDefined();
   const int *conn=_nodal_connec->getConstPointer();
   const int *connI=_nodal_connec_index->getConstPointer();
   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::size_t n=getAllGeoTypes().size();
+  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++)
     {
@@ -6240,7 +6601,7 @@ std::vector<int> MEDCouplingUMesh::getDistributionOfTypes() const throw(INTERP_K
  * If it exists a geometric type in \a this \b not in \a code \b no exception is thrown 
  * and a DataArrayInt instance is returned that the user has the responsability to deallocate.
  */
-DataArrayInt *MEDCouplingUMesh::checkTypeConsistencyAndContig(const std::vector<int>& code, const std::vector<const DataArrayInt *>& idsPerType) const throw(INTERP_KERNEL::Exception)
+DataArrayInt *MEDCouplingUMesh::checkTypeConsistencyAndContig(const std::vector<int>& code, const std::vector<const DataArrayInt *>& idsPerType) const
 {
   if(code.empty())
     throw INTERP_KERNEL::Exception("MEDCouplingUMesh::checkTypeConsistencyAndContig : code is empty, should not !");
@@ -6340,7 +6701,7 @@ DataArrayInt *MEDCouplingUMesh::checkTypeConsistencyAndContig(const std::vector<
  *              This vector can be empty in case of all geometric type cells are fully covered in ascending in the given input \a profile.
  * \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
  */
-void MEDCouplingUMesh::splitProfilePerType(const DataArrayInt *profile, std::vector<int>& code, std::vector<DataArrayInt *>& idsInPflPerType, std::vector<DataArrayInt *>& idsPerType) const throw(INTERP_KERNEL::Exception)
+void MEDCouplingUMesh::splitProfilePerType(const DataArrayInt *profile, std::vector<int>& code, std::vector<DataArrayInt *>& idsInPflPerType, std::vector<DataArrayInt *>& idsPerType) const
 {
   if(!profile)
     throw INTERP_KERNEL::Exception("MEDCouplingUMesh::splitProfilePerType : input profile is NULL !");
@@ -6417,7 +6778,7 @@ void MEDCouplingUMesh::splitProfilePerType(const DataArrayInt *profile, std::vec
  * The following equality should be verified 'nM1LevMesh->getMeshDimension()==this->getMeshDimension()-1'
  * This method returns 5+2 elements. 'desc', 'descIndx', 'revDesc', 'revDescIndx' and 'meshnM1' behaves exactly as ParaMEDMEM::MEDCouplingUMesh::buildDescendingConnectivity except the content as described after. The returned array specifies the n-1 mesh reordered by type as MEDMEM does. 'nM1LevMeshIds' contains the ids in returned 'meshnM1'. Finally 'meshnM1Old2New' contains numbering old2new that is to say the cell #k in coarse 'nM1LevMesh' will have the number ret[k] in returned mesh 'nM1LevMesh' MEDMEM reordered.
  */
-MEDCouplingUMesh *MEDCouplingUMesh::emulateMEDMEMBDC(const MEDCouplingUMesh *nM1LevMesh, DataArrayInt *desc, DataArrayInt *descIndx, DataArrayInt *&revDesc, DataArrayInt *&revDescIndx, DataArrayInt *& nM1LevMeshIds, DataArrayInt *&meshnM1Old2New) const throw(INTERP_KERNEL::Exception)
+MEDCouplingUMesh *MEDCouplingUMesh::emulateMEDMEMBDC(const MEDCouplingUMesh *nM1LevMesh, DataArrayInt *desc, DataArrayInt *descIndx, DataArrayInt *&revDesc, DataArrayInt *&revDescIndx, DataArrayInt *& nM1LevMeshIds, DataArrayInt *&meshnM1Old2New) const
 {
   checkFullyDefined();
   nM1LevMesh->checkFullyDefined();
@@ -6462,7 +6823,7 @@ MEDCouplingUMesh *MEDCouplingUMesh::emulateMEDMEMBDC(const MEDCouplingUMesh *nM1
  *          this array using decrRef() as it is no more needed.
  *  \throw If the nodal connectivity of cells is not defined.
  */
-DataArrayInt *MEDCouplingUMesh::sortCellsInMEDFileFrmt() throw(INTERP_KERNEL::Exception)
+DataArrayInt *MEDCouplingUMesh::sortCellsInMEDFileFrmt()
 {
   checkConnectivityFullyDefined();
   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret=getRenumArrForMEDFileFrmt();
@@ -6498,7 +6859,7 @@ bool MEDCouplingUMesh::checkConsecutiveCellTypes() const
  * 
  * \sa  MEDCouplingUMesh::checkConsecutiveCellTypesAndOrder
  */
-bool MEDCouplingUMesh::checkConsecutiveCellTypesForMEDFileFrmt() const throw(INTERP_KERNEL::Exception)
+bool MEDCouplingUMesh::checkConsecutiveCellTypesForMEDFileFrmt() const
 {
   return checkConsecutiveCellTypesAndOrder(MEDMEM_ORDER,MEDMEM_ORDER+N_MEDMEM_ORDER);
 }
@@ -6550,7 +6911,7 @@ bool MEDCouplingUMesh::checkConsecutiveCellTypesAndOrder(const INTERP_KERNEL::No
  * that tells for each cell the pos of its type in the array on type given in input parameter. The 2nd output parameter is an array with the same
  * number of tuples than input type array and with one component. This 2nd output array gives type by type the number of occurence of type in 'this'.
  */
-DataArrayInt *MEDCouplingUMesh::getLevArrPerCellTypes(const INTERP_KERNEL::NormalizedCellType *orderBg, const INTERP_KERNEL::NormalizedCellType *orderEnd, DataArrayInt *&nbPerType) const throw(INTERP_KERNEL::Exception)
+DataArrayInt *MEDCouplingUMesh::getLevArrPerCellTypes(const INTERP_KERNEL::NormalizedCellType *orderBg, const INTERP_KERNEL::NormalizedCellType *orderEnd, DataArrayInt *&nbPerType) const
 {
   checkConnectivityFullyDefined();
   int nbOfCells=getNumberOfCells();
@@ -6591,7 +6952,7 @@ DataArrayInt *MEDCouplingUMesh::getLevArrPerCellTypes(const INTERP_KERNEL::Norma
  *
  * \sa MEDCouplingUMesh::getRenumArrForConsecutiveCellTypesSpec, MEDCouplingUMesh::sortCellsInMEDFileFrmt.
  */
-DataArrayInt *MEDCouplingUMesh::getRenumArrForMEDFileFrmt() const throw(INTERP_KERNEL::Exception)
+DataArrayInt *MEDCouplingUMesh::getRenumArrForMEDFileFrmt() const
 {
   return getRenumArrForConsecutiveCellTypesSpec(MEDMEM_ORDER,MEDMEM_ORDER+N_MEDMEM_ORDER);
 }
@@ -6602,7 +6963,7 @@ DataArrayInt *MEDCouplingUMesh::getRenumArrForMEDFileFrmt() const throw(INTERP_K
  * The mesh after this call to MEDCouplingMesh::renumberCells will pass the test of MEDCouplingUMesh::checkConsecutiveCellTypesAndOrder with the same inputs.
  * The returned array minimizes the permutations that is to say the order of cells inside same geometric type remains the same.
  */
-DataArrayInt *MEDCouplingUMesh::getRenumArrForConsecutiveCellTypesSpec(const INTERP_KERNEL::NormalizedCellType *orderBg, const INTERP_KERNEL::NormalizedCellType *orderEnd) const throw(INTERP_KERNEL::Exception)
+DataArrayInt *MEDCouplingUMesh::getRenumArrForConsecutiveCellTypesSpec(const INTERP_KERNEL::NormalizedCellType *orderBg, const INTERP_KERNEL::NormalizedCellType *orderEnd) const
 {
   DataArrayInt *nbPerType=0;
   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> tmpa=getLevArrPerCellTypes(orderBg,orderEnd,nbPerType);
@@ -6655,7 +7016,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();
@@ -6687,42 +7048,99 @@ std::vector<MEDCouplingUMesh *> MEDCouplingUMesh::splitByType() const
  * \throw If the nodal connectivity of \a this is not fully defined.
  * \throw If the internal data is not coherent.
  */
-MEDCoupling1GTUMesh *MEDCouplingUMesh::convertIntoSingleGeoTypeMesh() const throw(INTERP_KERNEL::Exception)
+MEDCoupling1GTUMesh *MEDCouplingUMesh::convertIntoSingleGeoTypeMesh() const
 {
   checkConnectivityFullyDefined();
     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();
-  int typi=(int)typ;
   MEDCouplingAutoRefCountObjectPtr<MEDCoupling1GTUMesh> ret=MEDCoupling1GTUMesh::New(getName(),typ);
   ret->setCoords(getCoords());
   MEDCoupling1SGTUMesh *retC=dynamic_cast<MEDCoupling1SGTUMesh *>((MEDCoupling1GTUMesh*)ret);
   if(retC)
     {
-      int nbCells=getNumberOfCells();
-      int nbNodesPerCell=retC->getNumberOfNodesPerCell();
-      MEDCouplingAutoRefCountObjectPtr<DataArrayInt> connOut=DataArrayInt::New(); connOut->alloc(nbCells*nbNodesPerCell,1);
-      int *outPtr=connOut->getPointer();
-      const int *conn=_nodal_connec->begin();
-      const int *connI=_nodal_connec_index->begin();
-      nbNodesPerCell++;
-      for(int i=0;i<nbCells;i++,connI++)
-        {
-          if(conn[connI[0]]==typi && connI[1]-connI[0]==nbNodesPerCell)
-            outPtr=std::copy(conn+connI[0]+1,conn+connI[1],outPtr);
-          else
-            {
-              std::ostringstream oss; oss << "MEDCouplingUMesh::convertIntoSingleGeoTypeMesh : there something wrong in cell #" << i << " ! The type of cell is not those expected, or the length of nodal connectivity is not those expected (" << nbNodesPerCell-1 << ") !";
-              throw INTERP_KERNEL::Exception(oss.str().c_str());
-            }
-        }
-      retC->setNodalConnectivity(connOut);
+      MEDCouplingAutoRefCountObjectPtr<DataArrayInt> c=convertNodalConnectivityToStaticGeoTypeMesh();
+      retC->setNodalConnectivity(c);
     }
   else
-    throw INTERP_KERNEL::Exception("MEDCouplingUMesh::convertIntoSingleGeoTypeMesh : not implemented yet for non static geometric type !");
+    {
+      MEDCoupling1DGTUMesh *retD=dynamic_cast<MEDCoupling1DGTUMesh *>((MEDCoupling1GTUMesh*)ret);
+      if(!retD)
+        throw INTERP_KERNEL::Exception("MEDCouplingUMesh::convertIntoSingleGeoTypeMesh : Internal error !");
+      DataArrayInt *c=0,*ci=0;
+      convertNodalConnectivityToDynamicGeoTypeMesh(c,ci);
+      MEDCouplingAutoRefCountObjectPtr<DataArrayInt> cs(c),cis(ci);
+      retD->setNodalConnectivity(cs,cis);
+    }
   return ret.retn();
 }
 
+DataArrayInt *MEDCouplingUMesh::convertNodalConnectivityToStaticGeoTypeMesh() const
+{
+  checkConnectivityFullyDefined();
+    if(_types.size()!=1)
+    throw INTERP_KERNEL::Exception("MEDCouplingUMesh::convertNodalConnectivityToStaticGeoTypeMesh : current mesh does not contain exactly one geometric type !");
+  INTERP_KERNEL::NormalizedCellType typ=*_types.begin();
+  const INTERP_KERNEL::CellModel& cm=INTERP_KERNEL::CellModel::GetCellModel(typ);
+  if(cm.isDynamic())
+    {
+      std::ostringstream oss; oss << "MEDCouplingUMesh::convertNodalConnectivityToStaticGeoTypeMesh : this contains a single geo type (" << cm.getRepr() << ") but ";
+      oss << "this type is dynamic ! Only static geometric type is possible for that type ! call convertNodalConnectivityToDynamicGeoTypeMesh instead !";
+      throw INTERP_KERNEL::Exception(oss.str().c_str());
+    }
+  int nbCells=getNumberOfCells();
+  int typi=(int)typ;
+  int nbNodesPerCell=(int)cm.getNumberOfNodes();
+  MEDCouplingAutoRefCountObjectPtr<DataArrayInt> connOut=DataArrayInt::New(); connOut->alloc(nbCells*nbNodesPerCell,1);
+  int *outPtr=connOut->getPointer();
+  const int *conn=_nodal_connec->begin();
+  const int *connI=_nodal_connec_index->begin();
+  nbNodesPerCell++;
+  for(int i=0;i<nbCells;i++,connI++)
+    {
+      if(conn[connI[0]]==typi && connI[1]-connI[0]==nbNodesPerCell)
+        outPtr=std::copy(conn+connI[0]+1,conn+connI[1],outPtr);
+      else
+        {
+          std::ostringstream oss; oss << "MEDCouplingUMesh::convertNodalConnectivityToStaticGeoTypeMesh : there something wrong in cell #" << i << " ! The type of cell is not those expected, or the length of nodal connectivity is not those expected (" << nbNodesPerCell-1 << ") !";
+          throw INTERP_KERNEL::Exception(oss.str().c_str());
+        }
+    }
+  return connOut.retn();
+}
+
+void MEDCouplingUMesh::convertNodalConnectivityToDynamicGeoTypeMesh(DataArrayInt *&nodalConn, DataArrayInt *&nodalConnIndex) const
+{
+  static const char msg0[]="MEDCouplingUMesh::convertNodalConnectivityToDynamicGeoTypeMesh : nodal connectivity in this are invalid ! Call checkCoherency2 !";
+  checkConnectivityFullyDefined();
+  if(_types.size()!=1)
+    throw INTERP_KERNEL::Exception("MEDCouplingUMesh::convertNodalConnectivityToDynamicGeoTypeMesh : current mesh does not contain exactly one geometric type !");
+  int nbCells=getNumberOfCells(),lgth=_nodal_connec->getNumberOfTuples();
+  if(lgth<nbCells)
+    throw INTERP_KERNEL::Exception(msg0);
+  MEDCouplingAutoRefCountObjectPtr<DataArrayInt> c(DataArrayInt::New()),ci(DataArrayInt::New());
+  c->alloc(lgth-nbCells,1); ci->alloc(nbCells+1,1);
+  int *cp(c->getPointer()),*cip(ci->getPointer());
+  const int *incp(_nodal_connec->begin()),*incip(_nodal_connec_index->begin());
+  cip[0]=0;
+  for(int i=0;i<nbCells;i++,cip++,incip++)
+    {
+      int strt(incip[0]+1),stop(incip[1]);//+1 to skip geo type
+      int delta(stop-strt);
+      if(delta>=1)
+        {
+          if((strt>=0 && strt<lgth) && (stop>=0 && stop<=lgth))
+            cp=std::copy(incp+strt,incp+stop,cp);
+          else
+            throw INTERP_KERNEL::Exception(msg0);
+        }
+      else
+        throw INTERP_KERNEL::Exception(msg0);
+      cip[1]=cip[0]+delta;
+    }
+  nodalConn=c.retn(); nodalConnIndex=ci.retn();
+}
+
 /*!
  * This method takes in input a vector of MEDCouplingUMesh instances lying on the same coordinates with same mesh dimensions.
  * Each mesh in \b ms must be sorted by type with the same order (typically using MEDCouplingUMesh::sortCellsInMEDFileFrmt).
@@ -6794,7 +7212,7 @@ MEDCouplingUMesh *MEDCouplingUMesh::AggregateSortedByTypeMeshesOnSameCoords(cons
  * This method returns a newly created DataArrayInt instance.
  * This method retrieves cell ids in [ \a begin, \a end ) that have the type \a type.
  */
-DataArrayInt *MEDCouplingUMesh::keepCellIdsByType(INTERP_KERNEL::NormalizedCellType type, const int *begin, const int *end) const throw(INTERP_KERNEL::Exception)
+DataArrayInt *MEDCouplingUMesh::keepCellIdsByType(INTERP_KERNEL::NormalizedCellType type, const int *begin, const int *end) const
 {
   checkFullyDefined();
   const int *conn=_nodal_connec->getConstPointer();
@@ -6810,13 +7228,13 @@ DataArrayInt *MEDCouplingUMesh::keepCellIdsByType(INTERP_KERNEL::NormalizedCellT
  * This method makes the assumption that da->getNumberOfTuples()<this->getNumberOfCells(). This method makes the assumption that ids contained in 'da'
  * are in [0:getNumberOfCells())
  */
-DataArrayInt *MEDCouplingUMesh::convertCellArrayPerGeoType(const DataArrayInt *da) const throw(INTERP_KERNEL::Exception)
+DataArrayInt *MEDCouplingUMesh::convertCellArrayPerGeoType(const DataArrayInt *da) const
 {
   checkFullyDefined();
   const int *conn=_nodal_connec->getConstPointer();
   const int *connI=_nodal_connec_index->getConstPointer();
   int nbOfCells=getNumberOfCells();
-  std::set<INTERP_KERNEL::NormalizedCellType> types=getAllTypes();
+  std::set<INTERP_KERNEL::NormalizedCellType> types(getAllGeoTypes());
   int *tmp=new int[nbOfCells];
   for(std::set<INTERP_KERNEL::NormalizedCellType>::const_iterator iter=types.begin();iter!=types.end();iter++)
     {
@@ -6883,7 +7301,7 @@ MEDCouplingUMesh *MEDCouplingUMesh::keepSpecifiedCells(INTERP_KERNEL::Normalized
  * This method returns a vector of size 'this->getNumberOfCells()'.
  * This method retrieves for each cell in \a this if it is linear (false) or quadratic(true).
  */
-std::vector<bool> MEDCouplingUMesh::getQuadraticStatus() const throw(INTERP_KERNEL::Exception)
+std::vector<bool> MEDCouplingUMesh::getQuadraticStatus() const
 {
   int ncell=getNumberOfCells();
   std::vector<bool> ret(ncell);
@@ -6952,7 +7370,7 @@ DataArrayDouble *MEDCouplingUMesh::getBarycenterAndOwner() const
  * \throw If \a this is not fully defined (coordinates and connectivity)
  * \throw If there is presence in nodal connectivity in \a this of node ids not in [0, \c this->getNumberOfNodes() )
  */
-DataArrayDouble *MEDCouplingUMesh::computeIsoBarycenterOfNodesPerCell() const throw(INTERP_KERNEL::Exception)
+DataArrayDouble *MEDCouplingUMesh::computeIsoBarycenterOfNodesPerCell() const
 {
   checkFullyDefined();
   MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> ret=DataArrayDouble::New();
@@ -7029,8 +7447,10 @@ DataArrayDouble *MEDCouplingUMesh::computeIsoBarycenterOfNodesPerCell() const th
  *  \throw If the coordinates array is not set.
  *  \throw If the nodal connectivity of cells is not defined.
  *
+ *  \if ENABLE_EXAMPLES
  *  \ref cpp_mcumesh_getPartBarycenterAndOwner "Here is a C++ example".<br>
  *  \ref  py_mcumesh_getPartBarycenterAndOwner "Here is a Python example".
+ *  \endif
  */
 DataArrayDouble *MEDCouplingUMesh::getPartBarycenterAndOwner(const int *begin, const int *end) const
 {
@@ -7053,16 +7473,66 @@ DataArrayDouble *MEDCouplingUMesh::getPartBarycenterAndOwner(const int *begin, c
   return ret;
 }
 
+/*!
+ * Returns a DataArrayDouble instance giving for each cell in \a this the equation of plane given by "a*X+b*Y+c*Z+d=0".
+ * So the returned instance will have 4 components and \c this->getNumberOfCells() tuples.
+ * So this method expects that \a this has a spaceDimension equal to 3 and meshDimension equal to 2.
+ * The computation of the plane equation is done using each time the 3 first nodes of 2D cells.
+ * This method is useful to detect 2D cells in 3D space that are not coplanar.
+ * 
+ * \return DataArrayDouble * - a new instance of DataArrayDouble having 4 components and a number of tuples equal to number of cells in \a this.
+ * \throw If spaceDim!=3 or meshDim!=2.
+ * \throw If connectivity of \a this is invalid.
+ * \throw If connectivity of a cell in \a this points to an invalid node.
+ */
+DataArrayDouble *MEDCouplingUMesh::computePlaneEquationOf3DFaces() const
+{
+  MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> ret(DataArrayDouble::New());
+  int nbOfCells(getNumberOfCells()),nbOfNodes(getNumberOfNodes());
+  if(getSpaceDimension()!=3 || getMeshDimension()!=2)
+    throw INTERP_KERNEL::Exception("MEDCouplingUMesh::computePlaneEquationOf3DFaces : This method must be applied on a mesh having meshDimension equal 2 and a spaceDimension equal to 3 !");
+  ret->alloc(nbOfCells,4);
+  double *retPtr(ret->getPointer());
+  const int *nodal(_nodal_connec->begin()),*nodalI(_nodal_connec_index->begin());
+  const double *coor(_coords->begin());
+  for(int i=0;i<nbOfCells;i++,nodalI++,retPtr+=4)
+    {
+      double matrix[16]={0,0,0,1,0,0,0,1,0,0,0,1,1,1,1,0},matrix2[16];
+      if(nodalI[1]-nodalI[0]>=3)
+        {
+          for(int j=0;j<3;j++)
+            {
+              int nodeId(nodal[nodalI[0]+1+j]);
+              if(nodeId>=0 && nodeId<nbOfNodes)
+                std::copy(coor+nodeId*3,coor+(nodeId+1)*3,matrix+4*j);
+              else
+                {
+                  std::ostringstream oss; oss << "MEDCouplingUMesh::computePlaneEquationOf3DFaces : invalid 2D cell #" << i << " ! This cell points to an invalid nodeId : " << nodeId << " !";
+                  throw INTERP_KERNEL::Exception(oss.str().c_str());
+                }
+            }
+        }
+      else
+        {
+          std::ostringstream oss; oss << "MEDCouplingUMesh::computePlaneEquationOf3DFaces : invalid 2D cell #" << i << " ! Must be constitued by more than 3 nodes !";
+          throw INTERP_KERNEL::Exception(oss.str().c_str());
+        }
+      INTERP_KERNEL::inverseMatrix(matrix,4,matrix2);
+      retPtr[0]=matrix2[3]; retPtr[1]=matrix2[7]; retPtr[2]=matrix2[11]; retPtr[3]=matrix2[15];
+    }
+  return ret.retn();
+}
+
 /*!
  * This method expects as input a DataArrayDouble non nul instance 'da' that should be allocated. If not an exception is thrown.
  * 
  */
-MEDCouplingUMesh *MEDCouplingUMesh::Build0DMeshFromCoords(DataArrayDouble *da) throw(INTERP_KERNEL::Exception)
+MEDCouplingUMesh *MEDCouplingUMesh::Build0DMeshFromCoords(DataArrayDouble *da)
 {
   if(!da)
     throw INTERP_KERNEL::Exception("MEDCouplingUMesh::Build0DMeshFromCoords : instance of DataArrayDouble must be not null !");
   da->checkAllocated();
-  MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> ret=MEDCouplingUMesh::New(da->getName().c_str(),0);
+  MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> ret=MEDCouplingUMesh::New(da->getName(),0);
   ret->setCoords(da);
   int nbOfTuples=da->getNumberOfTuples();
   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> c=DataArrayInt::New();
@@ -7095,7 +7565,7 @@ MEDCouplingUMesh *MEDCouplingUMesh::Build0DMeshFromCoords(DataArrayDouble *da) t
  *  \throw If \a mesh1->getMeshDimension() < 0 or \a mesh2->getMeshDimension() < 0.
  *  \throw If \a mesh1->getMeshDimension() != \a mesh2->getMeshDimension().
  */
-MEDCouplingUMesh *MEDCouplingUMesh::MergeUMeshes(const MEDCouplingUMesh *mesh1, const MEDCouplingUMesh *mesh2) throw(INTERP_KERNEL::Exception)
+MEDCouplingUMesh *MEDCouplingUMesh::MergeUMeshes(const MEDCouplingUMesh *mesh1, const MEDCouplingUMesh *mesh2)
 {
   std::vector<const MEDCouplingUMesh *> tmp(2);
   tmp[0]=const_cast<MEDCouplingUMesh *>(mesh1); tmp[1]=const_cast<MEDCouplingUMesh *>(mesh2);
@@ -7116,7 +7586,7 @@ MEDCouplingUMesh *MEDCouplingUMesh::MergeUMeshes(const MEDCouplingUMesh *mesh1,
  *  \throw If \a a[ *i* ]->getMeshDimension() < 0.
  *  \throw If the meshes in \a a are of different dimension (getMeshDimension()).
 */
-MEDCouplingUMesh *MEDCouplingUMesh::MergeUMeshes(std::vector<const MEDCouplingUMesh *>& a) throw(INTERP_KERNEL::Exception)
+MEDCouplingUMesh *MEDCouplingUMesh::MergeUMeshes(std::vector<const MEDCouplingUMesh *>& a)
 {
   std::size_t sz=a.size();
   if(sz==0)
@@ -7149,7 +7619,7 @@ MEDCouplingUMesh *MEDCouplingUMesh::MergeUMeshes(std::vector<const MEDCouplingUM
 
 /// @cond INTERNAL
 
-MEDCouplingUMesh *MEDCouplingUMesh::MergeUMeshesLL(std::vector<const MEDCouplingUMesh *>& a) throw(INTERP_KERNEL::Exception)
+MEDCouplingUMesh *MEDCouplingUMesh::MergeUMeshesLL(std::vector<const MEDCouplingUMesh *>& a)
 {
   if(a.empty())
     throw INTERP_KERNEL::Exception("MEDCouplingUMesh::MergeUMeshes : input array must be NON EMPTY !");
@@ -7221,7 +7691,7 @@ MEDCouplingUMesh *MEDCouplingUMesh::MergeUMeshesLL(std::vector<const MEDCoupling
  *  \throw If \a mesh1->getMeshDimension() < 0 or \a mesh2->getMeshDimension() < 0.
  *  \throw If \a mesh1->getMeshDimension() != \a mesh2->getMeshDimension().
  */
-MEDCouplingUMesh *MEDCouplingUMesh::MergeUMeshesOnSameCoords(const MEDCouplingUMesh *mesh1, const MEDCouplingUMesh *mesh2) throw(INTERP_KERNEL::Exception)
+MEDCouplingUMesh *MEDCouplingUMesh::MergeUMeshesOnSameCoords(const MEDCouplingUMesh *mesh1, const MEDCouplingUMesh *mesh2)
 {
   std::vector<const MEDCouplingUMesh *> tmp(2);
   tmp[0]=mesh1; tmp[1]=mesh2;
@@ -7354,7 +7824,7 @@ MEDCouplingUMesh *MEDCouplingUMesh::FuseUMeshesOnSameCoords(const std::vector<co
  *  \throw If the nodal connectivity of cells is not defined in any of \a meshes.
  *  \throw If \a meshes are of different space dimension.
  */
-void MEDCouplingUMesh::PutUMeshesOnSameAggregatedCoords(const std::vector<MEDCouplingUMesh *>& meshes) throw(INTERP_KERNEL::Exception)
+void MEDCouplingUMesh::PutUMeshesOnSameAggregatedCoords(const std::vector<MEDCouplingUMesh *>& meshes)
 {
   std::size_t sz=meshes.size();
   if(sz==0 || sz==1)
@@ -7408,7 +7878,7 @@ void MEDCouplingUMesh::PutUMeshesOnSameAggregatedCoords(const std::vector<MEDCou
  *  \throw If the \a meshes do not share the same node coordinates array.
  *  \throw If the nodal connectivity of cells is not defined in any of \a meshes.
  */
-void MEDCouplingUMesh::MergeNodesOnUMeshesSharingSameCoords(const std::vector<MEDCouplingUMesh *>& meshes, double eps) throw(INTERP_KERNEL::Exception)
+void MEDCouplingUMesh::MergeNodesOnUMeshesSharingSameCoords(const std::vector<MEDCouplingUMesh *>& meshes, double eps)
 {
   if(meshes.empty())
     return ;
@@ -7536,17 +8006,43 @@ void MEDCouplingUMesh::AppendExtrudedCell(const int *connBg, const int *connEnd,
  */
 bool MEDCouplingUMesh::IsPolygonWellOriented(bool isQuadratic, const double *vec, const int *begin, const int *end, const double *coords)
 {
+  std::size_t i, ip1;
   double v[3]={0.,0.,0.};
   std::size_t sz=std::distance(begin,end);
   if(isQuadratic)
     sz/=2;
-  for(std::size_t i=0;i<sz;i++)
+  for(i=0;i<sz;i++)
     {
       v[0]+=coords[3*begin[i]+1]*coords[3*begin[(i+1)%sz]+2]-coords[3*begin[i]+2]*coords[3*begin[(i+1)%sz]+1];
       v[1]+=coords[3*begin[i]+2]*coords[3*begin[(i+1)%sz]]-coords[3*begin[i]]*coords[3*begin[(i+1)%sz]+2];
       v[2]+=coords[3*begin[i]]*coords[3*begin[(i+1)%sz]+1]-coords[3*begin[i]+1]*coords[3*begin[(i+1)%sz]];
     }
-  return vec[0]*v[0]+vec[1]*v[1]+vec[2]*v[2]>0.;
+  double ret = vec[0]*v[0]+vec[1]*v[1]+vec[2]*v[2];
+
+  // Try using quadratic points if standard points are degenerated (for example a QPOLYG with two
+  // SEG3 forming a circle):
+  if (fabs(ret) < INTERP_KERNEL::DEFAULT_ABS_TOL && isQuadratic)
+    {
+      v[0] = 0.0; v[1] = 0.0; v[2] = 0.0;
+      for(std::size_t j=0;j<sz;j++)
+        {
+          if (j%2)  // current point i is quadratic, next point i+1 is standard
+            {
+              i = sz+j;
+              ip1 = (j+1)%sz; // ip1 = "i+1"
+            }
+          else      // current point i is standard, next point i+1 is quadratic
+            {
+              i = j;
+              ip1 = j+sz;
+            }
+          v[0]+=coords[3*begin[i]+1]*coords[3*begin[ip1]+2]-coords[3*begin[i]+2]*coords[3*begin[ip1]+1];
+          v[1]+=coords[3*begin[i]+2]*coords[3*begin[ip1]]-coords[3*begin[i]]*coords[3*begin[ip1]+2];
+          v[2]+=coords[3*begin[i]]*coords[3*begin[ip1]+1]-coords[3*begin[i]+1]*coords[3*begin[ip1]];
+        }
+      ret = vec[0]*v[0]+vec[1]*v[1]+vec[2]*v[2];
+    }
+  return (ret>0.);
 }
 
 /*!
@@ -7636,7 +8132,7 @@ bool MEDCouplingUMesh::IsPyra5WellOriented(const int *begin, const int *end, con
  * \param [in] end end of nodal connectivity of a single polyhedron cell (excluded)
  * \param [out] res the result is put at the end of the vector without any alteration of the data.
  */
-void MEDCouplingUMesh::SimplifyPolyhedronCell(double eps, const DataArrayDouble *coords, const int *begin, const int *end, DataArrayInt *res) throw(INTERP_KERNEL::Exception)
+void MEDCouplingUMesh::SimplifyPolyhedronCell(double eps, const DataArrayDouble *coords, const int *begin, const int *end, DataArrayInt *res)
 {
   int nbFaces=std::count(begin+1,end,-1)+1;
   MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> v=DataArrayDouble::New(); v->alloc(nbFaces,3);
@@ -7728,7 +8224,7 @@ void MEDCouplingUMesh::SimplifyPolyhedronCell(double eps, const DataArrayDouble
  * \param [out] v the normalized vector of size 3
  * \param [out] p the pos of plane
  */
-void MEDCouplingUMesh::ComputeVecAndPtOfFace(double eps, const double *coords, const int *begin, const int *end, double *v, double *p) throw(INTERP_KERNEL::Exception)
+void MEDCouplingUMesh::ComputeVecAndPtOfFace(double eps, const double *coords, const int *begin, const int *end, double *v, double *p)
 {
   std::size_t nbPoints=std::distance(begin,end);
   if(nbPoints<3)
@@ -7774,7 +8270,7 @@ void MEDCouplingUMesh::ComputeVecAndPtOfFace(double eps, const double *coords, c
  * This method tries to obtain a well oriented polyhedron.
  * If the algorithm fails, an exception will be thrown.
  */
-void MEDCouplingUMesh::TryToCorrectPolyhedronOrientation(int *begin, int *end, const double *coords) throw(INTERP_KERNEL::Exception)
+void MEDCouplingUMesh::TryToCorrectPolyhedronOrientation(int *begin, int *end, const double *coords)
 {
   std::list< std::pair<int,int> > edgesOK,edgesFinished;
   std::size_t nbOfFaces=std::count(begin,end,-1)+1;
@@ -7845,34 +8341,20 @@ void MEDCouplingUMesh::TryToCorrectPolyhedronOrientation(int *begin, int *end, c
     }
 }
 
-/*!
- * This method makes the assumption spacedimension == meshdimension == 2.
- * This method works only for linear cells.
- * 
- * \return a newly allocated array containing the connectivity of a polygon type enum included (NORM_POLYGON in pos#0)
- */
-DataArrayInt *MEDCouplingUMesh::buildUnionOf2DMesh() const throw(INTERP_KERNEL::Exception)
+DataArrayInt *MEDCouplingUMesh::buildUnionOf2DMeshLinear(const MEDCouplingUMesh *skin, const DataArrayInt *n2o) const
 {
-  if(getMeshDimension()!=2 || getSpaceDimension()!=2)
-    throw INTERP_KERNEL::Exception("MEDCouplingUMesh::buildUnionOf2DMesh : meshdimension, spacedimension must be equal to 2 !");
-  MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> m=computeSkin();
-  MEDCouplingAutoRefCountObjectPtr<DataArrayInt> o2n=m->zipCoordsTraducer();
-  int nbOfNodesExpected=m->getNumberOfNodes();
-  if(m->getNumberOfCells()!=nbOfNodesExpected)
-    throw INTERP_KERNEL::Exception("MEDCouplingUMesh::buildUnionOf2DMesh : the mesh 2D in input appears to be not in a single part or a quadratic 2D mesh !");
-  MEDCouplingAutoRefCountObjectPtr<DataArrayInt> n2o=o2n->invertArrayO2N2N2O(m->getNumberOfNodes());
-  const int *n2oPtr=n2o->getConstPointer();
+  int nbOfNodesExpected(skin->getNumberOfNodes());
+  const int *n2oPtr(n2o->getConstPointer());
   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> revNodal(DataArrayInt::New()),revNodalI(DataArrayInt::New());
-  m->getReverseNodalConnectivity(revNodal,revNodalI);
-  const int *revNodalPtr=revNodal->getConstPointer(),*revNodalIPtr=revNodalI->getConstPointer();
-  const int *nodalPtr=m->getNodalConnectivity()->getConstPointer();
-  const int *nodalIPtr=m->getNodalConnectivityIndex()->getConstPointer();
-  MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret=DataArrayInt::New(); ret->alloc(nbOfNodesExpected+1,1);
-  int *work=ret->getPointer();  *work++=INTERP_KERNEL::NORM_POLYGON;
+  skin->getReverseNodalConnectivity(revNodal,revNodalI);
+  const int *revNodalPtr(revNodal->getConstPointer()),*revNodalIPtr(revNodalI->getConstPointer());
+  const int *nodalPtr(skin->getNodalConnectivity()->getConstPointer());
+  const int *nodalIPtr(skin->getNodalConnectivityIndex()->getConstPointer());
+  MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret(DataArrayInt::New()); ret->alloc(nbOfNodesExpected+1,1);
+  int *work(ret->getPointer());  *work++=INTERP_KERNEL::NORM_POLYGON;
   if(nbOfNodesExpected<1)
     return ret.retn();
-  int prevCell=0;
-  int prevNode=nodalPtr[nodalIPtr[0]+1];
+  int prevCell(0),prevNode(nodalPtr[nodalIPtr[0]+1]);
   *work++=n2oPtr[prevNode];
   for(int i=1;i<nbOfNodesExpected;i++)
     {
@@ -7882,7 +8364,7 @@ DataArrayInt *MEDCouplingUMesh::buildUnionOf2DMesh() const throw(INTERP_KERNEL::
           conn.erase(prevNode);
           if(conn.size()==1)
             {
-              int curNode=*(conn.begin());
+              int curNode(*(conn.begin()));
               *work++=n2oPtr[curNode];
               std::set<int> shar(revNodalPtr+revNodalIPtr[curNode],revNodalPtr+revNodalIPtr[curNode+1]);
               shar.erase(prevCell);
@@ -7892,24 +8374,96 @@ DataArrayInt *MEDCouplingUMesh::buildUnionOf2DMesh() const throw(INTERP_KERNEL::
                   prevNode=curNode;
                 }
               else
-                throw INTERP_KERNEL::Exception("MEDCouplingUMesh::buildUnionOf2DMesh : presence of unexpected 2 !");
+                throw INTERP_KERNEL::Exception("MEDCouplingUMesh::buildUnionOf2DMeshLinear : presence of unexpected 2 !");
             }
           else
-            throw INTERP_KERNEL::Exception("MEDCouplingUMesh::buildUnionOf2DMesh : presence of unexpected 1 !");
+            throw INTERP_KERNEL::Exception("MEDCouplingUMesh::buildUnionOf2DMeshLinear : presence of unexpected 1 !");
         }
       else
-        throw INTERP_KERNEL::Exception("MEDCouplingUMesh::buildUnionOf2DMesh : presence of unexpected cell !");
+        throw INTERP_KERNEL::Exception("MEDCouplingUMesh::buildUnionOf2DMeshLinear : presence of unexpected cell !");
     }
   return ret.retn();
 }
 
+DataArrayInt *MEDCouplingUMesh::buildUnionOf2DMeshQuadratic(const MEDCouplingUMesh *skin, const DataArrayInt *n2o) const
+{
+  int nbOfNodesExpected(skin->getNumberOfNodes());
+  int nbOfTurn(nbOfNodesExpected/2);
+  const int *n2oPtr(n2o->getConstPointer());
+  MEDCouplingAutoRefCountObjectPtr<DataArrayInt> revNodal(DataArrayInt::New()),revNodalI(DataArrayInt::New());
+  skin->getReverseNodalConnectivity(revNodal,revNodalI);
+  const int *revNodalPtr(revNodal->getConstPointer()),*revNodalIPtr(revNodalI->getConstPointer());
+  const int *nodalPtr(skin->getNodalConnectivity()->getConstPointer());
+  const int *nodalIPtr(skin->getNodalConnectivityIndex()->getConstPointer());
+  MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret(DataArrayInt::New()); ret->alloc(nbOfNodesExpected+1,1);
+  int *work(ret->getPointer());  *work++=INTERP_KERNEL::NORM_QPOLYG;
+  if(nbOfNodesExpected<1)
+    return ret.retn();
+  int prevCell(0),prevNode(nodalPtr[nodalIPtr[0]+1]);
+  *work=n2oPtr[prevNode]; work[nbOfTurn]=n2oPtr[nodalPtr[nodalIPtr[0]+3]]; work++;
+  for(int i=1;i<nbOfTurn;i++)
+    {
+      if(nodalIPtr[prevCell+1]-nodalIPtr[prevCell]==4)
+        {
+          std::set<int> conn(nodalPtr+nodalIPtr[prevCell]+1,nodalPtr+nodalIPtr[prevCell]+3);
+          conn.erase(prevNode);
+          if(conn.size()==1)
+            {
+              int curNode(*(conn.begin()));
+              *work=n2oPtr[curNode];
+              std::set<int> shar(revNodalPtr+revNodalIPtr[curNode],revNodalPtr+revNodalIPtr[curNode+1]);
+              shar.erase(prevCell);
+              if(shar.size()==1)
+                {
+                  int curCell(*(shar.begin()));
+                  work[nbOfTurn]=n2oPtr[nodalPtr[nodalIPtr[curCell]+3]];
+                  prevCell=curCell;
+                  prevNode=curNode;
+                  work++;
+                }
+              else
+                throw INTERP_KERNEL::Exception("MEDCouplingUMesh::buildUnionOf2DMeshQuadratic : presence of unexpected 2 !");
+            }
+          else
+            throw INTERP_KERNEL::Exception("MEDCouplingUMesh::buildUnionOf2DMeshQuadratic : presence of unexpected 1 !");
+        }
+      else
+        throw INTERP_KERNEL::Exception("MEDCouplingUMesh::buildUnionOf2DMeshQuadratic : presence of unexpected cell !");
+    }
+  return ret.retn();
+}
+
+/*!
+ * This method makes the assumption spacedimension == meshdimension == 2.
+ * This method works only for linear cells.
+ * 
+ * \return a newly allocated array containing the connectivity of a polygon type enum included (NORM_POLYGON in pos#0)
+ */
+DataArrayInt *MEDCouplingUMesh::buildUnionOf2DMesh() const
+{
+  if(getMeshDimension()!=2 || getSpaceDimension()!=2)
+    throw INTERP_KERNEL::Exception("MEDCouplingUMesh::buildUnionOf2DMesh : meshdimension, spacedimension must be equal to 2 !");
+  MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> skin(computeSkin());
+  int oldNbOfNodes(skin->getNumberOfNodes());
+  MEDCouplingAutoRefCountObjectPtr<DataArrayInt> o2n(skin->zipCoordsTraducer());
+  int nbOfNodesExpected(skin->getNumberOfNodes());
+  MEDCouplingAutoRefCountObjectPtr<DataArrayInt> n2o(o2n->invertArrayO2N2N2O(oldNbOfNodes));
+  int nbCells(skin->getNumberOfCells());
+  if(nbCells==nbOfNodesExpected)
+    return buildUnionOf2DMeshLinear(skin,n2o);
+  else if(2*nbCells==nbOfNodesExpected)
+    return buildUnionOf2DMeshQuadratic(skin,n2o);
+  else
+    throw INTERP_KERNEL::Exception("MEDCouplingUMesh::buildUnionOf2DMesh : the mesh 2D in input appears to be not in a single part of a 2D mesh !");
+}
+
 /*!
  * This method makes the assumption spacedimension == meshdimension == 3.
  * This method works only for linear cells.
  * 
  * \return a newly allocated array containing the connectivity of a polygon type enum included (NORM_POLYHED in pos#0)
  */
-DataArrayInt *MEDCouplingUMesh::buildUnionOf3DMesh() const throw(INTERP_KERNEL::Exception)
+DataArrayInt *MEDCouplingUMesh::buildUnionOf3DMesh() const
 {
   if(getMeshDimension()!=3 || getSpaceDimension()!=3)
     throw INTERP_KERNEL::Exception("MEDCouplingUMesh::buildUnionOf3DMesh : meshdimension, spacedimension must be equal to 2 !");
@@ -7934,7 +8488,7 @@ DataArrayInt *MEDCouplingUMesh::buildUnionOf3DMesh() const throw(INTERP_KERNEL::
  * This method put in zip format into parameter 'zipFrmt' in full interlace mode.
  * This format is often asked by INTERP_KERNEL algorithms to avoid many indirections into coordinates array.
  */
-void MEDCouplingUMesh::FillInCompact3DMode(int spaceDim, int nbOfNodesInCell, const int *conn, const double *coo, double *zipFrmt) throw(INTERP_KERNEL::Exception)
+void MEDCouplingUMesh::FillInCompact3DMode(int spaceDim, int nbOfNodesInCell, const int *conn, const double *coo, double *zipFrmt)
 {
   double *w=zipFrmt;
   if(spaceDim==3)
@@ -7952,12 +8506,12 @@ 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
 {
   int nbOfCells=getNumberOfCells();
   if(nbOfCells<=0)
     throw INTERP_KERNEL::Exception("MEDCouplingUMesh::writeVTK : the unstructured mesh has no cells !");
-  static const int PARAMEDMEM2VTKTYPETRADUCER[INTERP_KERNEL::NORM_MAXTYPE+1]={1,3,21,5,9,7,22,34,23,28,-1,-1,-1,-1,10,14,13,-1,12,-1,24,-1,16,27,-1,26,-1,29,-1,-1,25,42,-1,4};
+  static const int PARAMEDMEM2VTKTYPETRADUCER[INTERP_KERNEL::NORM_MAXTYPE+1]={1,3,21,5,9,7,22,34,23,28,-1,-1,-1,-1,10,14,13,-1,12,-1,24,-1,16,27,-1,26,-1,29,-1,-1,25,42,36,4};
   ofs << "  <" << getVTKDataSetType() << ">\n";
   ofs << "    <Piece NumberOfPoints=\"" << getNumberOfNodes() << "\" NumberOfCells=\"" << nbOfCells << "\">\n";
   ofs << "      <PointData>\n" << pointData << std::endl;
@@ -7966,11 +8520,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";
@@ -8001,12 +8555,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++)
@@ -8023,15 +8577,15 @@ 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";
 }
 
-void MEDCouplingUMesh::reprQuickOverview(std::ostream& stream) const throw(INTERP_KERNEL::Exception)
+void MEDCouplingUMesh::reprQuickOverview(std::ostream& stream) const
 {
   stream << "MEDCouplingUMesh C++ instance at " << this << ". Name : \"" << getName() << "\".";
   if(_mesh_dim==-2)
@@ -8056,14 +8610,17 @@ void MEDCouplingUMesh::reprQuickOverview(std::ostream& stream) const throw(INTER
   stream << std::endl << "Number of cells : " << lgth-1 << ".";
 }
 
-std::string MEDCouplingUMesh::getVTKDataSetType() const throw(INTERP_KERNEL::Exception)
+std::string MEDCouplingUMesh::getVTKDataSetType() const
 {
   return std::string("UnstructuredGrid");
 }
 
 /*!
  * Partitions the first given 2D mesh using the second given 2D mesh as a tool, and
- * returns a result mesh constituted by polygons. The meshes should be in 2D space. In
+ * returns a result mesh constituted by polygons.
+ * Thus the final result contains all nodes from m1 plus new nodes. However it doesn't necessarily contains
+ * all nodes from m2.
+ * The meshes should be in 2D space. In
  * addition, returns two arrays mapping cells of the result mesh to cells of the input
  * meshes.
  *  \param [in] m1 - the first input mesh which is a partitioned object.
@@ -8084,31 +8641,40 @@ std::string MEDCouplingUMesh::getVTKDataSetType() const throw(INTERP_KERNEL::Exc
  *  \throw If the nodal connectivity of cells is not defined in any of the meshes.
  *  \throw If any of the meshes is not a 2D mesh in 2D space.
  */
-MEDCouplingUMesh *MEDCouplingUMesh::Intersect2DMeshes(const MEDCouplingUMesh *m1, const MEDCouplingUMesh *m2, double eps, DataArrayInt *&cellNb1, DataArrayInt *&cellNb2) throw(INTERP_KERNEL::Exception)
+MEDCouplingUMesh *MEDCouplingUMesh::Intersect2DMeshes(const MEDCouplingUMesh *m1, const MEDCouplingUMesh *m2,
+                                                      double eps, DataArrayInt *&cellNb1, DataArrayInt *&cellNb2)
 {
   m1->checkFullyDefined();
   m2->checkFullyDefined();
   if(m1->getMeshDimension()!=2 || m1->getSpaceDimension()!=2 || m2->getMeshDimension()!=2 || m2->getSpaceDimension()!=2)
     throw INTERP_KERNEL::Exception("MEDCouplingUMesh::Intersect2DMeshes works on umeshes m1 AND m2  with meshdim equal to 2 and spaceDim equal to 2 too!");
+
+  // Step 1: compute all edge intersections (new nodes)
   std::vector< std::vector<int> > intersectEdge1, colinear2, subDiv2;
-  MEDCouplingUMesh *m1Desc=0,*m2Desc=0;
+  MEDCouplingUMesh *m1Desc=0,*m2Desc=0; // descending connec. meshes
   DataArrayInt *desc1=0,*descIndx1=0,*revDesc1=0,*revDescIndx1=0,*desc2=0,*descIndx2=0,*revDesc2=0,*revDescIndx2=0;
-  std::vector<double> addCoo,addCoordsQuadratic;
+  std::vector<double> addCoo,addCoordsQuadratic;  // coordinates of newly created nodes
   INTERP_KERNEL::QUADRATIC_PLANAR::_precision=eps;
   INTERP_KERNEL::QUADRATIC_PLANAR::_arc_detection_precision=eps;
-  IntersectDescending2DMeshes(m1,m2,eps,intersectEdge1,colinear2, subDiv2,m1Desc,desc1,descIndx1,revDesc1,revDescIndx1,
-                              m2Desc,desc2,descIndx2,revDesc2,revDescIndx2,addCoo);
+  IntersectDescending2DMeshes(m1,m2,eps,intersectEdge1,colinear2, subDiv2,
+                                      m1Desc,desc1,descIndx1,revDesc1,revDescIndx1,
+                              addCoo, m2Desc,desc2,descIndx2,revDesc2,revDescIndx2);
   revDesc1->decrRef(); revDescIndx1->decrRef(); revDesc2->decrRef(); revDescIndx2->decrRef();
   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> dd1(desc1),dd2(descIndx1),dd3(desc2),dd4(descIndx2);
   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> dd5(m1Desc),dd6(m2Desc);
+
+  // Step 2: re-order newly created nodes according to the ordering found in m2
   std::vector< std::vector<int> > intersectEdge2;
   BuildIntersectEdges(m1Desc,m2Desc,addCoo,subDiv2,intersectEdge2);
   subDiv2.clear(); dd5=0; dd6=0;
+
+  // Step 3:
   std::vector<int> cr,crI; //no DataArrayInt because interface with Geometric2D
   std::vector<int> cNb1,cNb2; //no DataArrayInt because interface with Geometric2D
   BuildIntersecting2DCellsFromEdges(eps,m1,desc1->getConstPointer(),descIndx1->getConstPointer(),intersectEdge1,colinear2,m2,desc2->getConstPointer(),descIndx2->getConstPointer(),intersectEdge2,addCoo,
                                     /* outputs -> */addCoordsQuadratic,cr,crI,cNb1,cNb2);
-  //
+
+  // Step 4: Prepare final result:
   MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> addCooDa=DataArrayDouble::New();
   addCooDa->alloc((int)(addCoo.size())/2,2);
   std::copy(addCoo.begin(),addCoo.end(),addCooDa->getPointer());
@@ -8129,6 +8695,15 @@ MEDCouplingUMesh *MEDCouplingUMesh::Intersect2DMeshes(const MEDCouplingUMesh *m1
   return ret.retn();
 }
 
+
+/**
+ * Private. Third step of the partitioning algorithm (Intersect2DMeshes): reconstruct full 2D cells from the
+ * (newly created) nodes corresponding to the edge intersections.
+ * Output params:
+ * @param[out] cr, crI connectivity of the resulting mesh
+ * @param[out] cNb1, cNb2 correspondance arrays giving for the merged mesh the initial cells IDs in m1 / m2
+ * TODO: describe input parameters
+ */
 void MEDCouplingUMesh::BuildIntersecting2DCellsFromEdges(double eps, const MEDCouplingUMesh *m1, const int *desc1, const int *descIndx1,
                                                          const std::vector<std::vector<int> >& intesctEdges1, const std::vector< std::vector<int> >& colinear2,
                                                          const MEDCouplingUMesh *m2, const int *desc2, const int *descIndx2, const std::vector<std::vector<int> >& intesctEdges2,
@@ -8136,7 +8711,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();
@@ -8146,21 +8720,24 @@ 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());
+  // Here a BBTree on 2D-cells, not on segments:
+  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;
       INTERP_KERNEL::NormalizedCellType typ=(INTERP_KERNEL::NormalizedCellType)conn1[connI1[i]];
       const INTERP_KERNEL::CellModel& cm=INTERP_KERNEL::CellModel::GetCellModel(typ);
+      // Populate mapp and mappRev with nodes from the current cell (i) from mesh1 - this also builds the Node* objects:
       MEDCouplingUMeshBuildQPFromMesh3(coo1,offset1,coo2,offset2,addCoords,desc1+descIndx1[i],desc1+descIndx1[i+1],intesctEdges1,/* output */mapp,mappRev);
+      // pol1 is the full cell from mesh2, in QP format, with all the additional intersecting nodes.
       pol1.buildFromCrudeDataArray(mappRev,cm.isQuadratic(),conn1+connI1[i]+1,coo1,
                                    desc1+descIndx1[i],desc1+descIndx1[i+1],intesctEdges1);
       //
@@ -8170,16 +8747,18 @@ void MEDCouplingUMesh::BuildIntersecting2DCellsFromEdges(double eps, const MEDCo
       for(it1.first();!it1.finished();it1.next())
         edges1.insert(it1.current()->getPtr());
       //
-      std::map<int,std::vector<INTERP_KERNEL::ElementaryEdge *> > edgesIn2ForShare;
+      std::map<int,std::vector<INTERP_KERNEL::ElementaryEdge *> > edgesIn2ForShare; // common edges
       std::vector<INTERP_KERNEL::QuadraticPolygon> pol2s(candidates2.size());
       int ii=0;
       for(std::vector<int>::const_iterator it2=candidates2.begin();it2!=candidates2.end();it2++,ii++)
         {
           INTERP_KERNEL::NormalizedCellType typ2=(INTERP_KERNEL::NormalizedCellType)conn2[connI2[*it2]];
           const INTERP_KERNEL::CellModel& cm2=INTERP_KERNEL::CellModel::GetCellModel(typ2);
+          // Complete mapping with elements coming from the current cell it2 in mesh2:
           MEDCouplingUMeshBuildQPFromMesh3(coo1,offset1,coo2,offset2,addCoords,desc2+descIndx2[*it2],desc2+descIndx2[*it2+1],intesctEdges2,/* output */mapp,mappRev);
+          // pol2 is the new QP in the final merged result.
           pol2s[ii].buildFromCrudeDataArray2(mappRev,cm2.isQuadratic(),conn2+connI2[*it2]+1,coo2,desc2+descIndx2[*it2],desc2+descIndx2[*it2+1],intesctEdges2,
-                                             pol1,desc1+descIndx1[i],desc1+descIndx1[i+1],intesctEdges1,colinear2,edgesIn2ForShare);
+                                             pol1,desc1+descIndx1[i],desc1+descIndx1[i+1],intesctEdges1,colinear2, /* output */ edgesIn2ForShare);
         }
       ii=0;
       for(std::vector<int>::const_iterator it2=candidates2.begin();it2!=candidates2.end();it2++,ii++)
@@ -8189,6 +8768,8 @@ void MEDCouplingUMesh::BuildIntersecting2DCellsFromEdges(double eps, const MEDCo
           //MEDCouplingUMeshAssignOnLoc(pol1,pol2,desc1+descIndx1[i],desc1+descIndx1[i+1],intesctEdges1,desc2+descIndx2[*it2],desc2+descIndx2[*it2+1],intesctEdges2,colinear2);
           pol1.buildPartitionsAbs(pol2s[ii],edges1,edgesBoundary2,mapp,i,*it2,offset3,addCoordsQuadratic,cr,crI,cNb1,cNb2);
         }
+      // Deals with remaining (non-consumed) edges from m1: these are the edges that were never touched
+      // by m2 but that we still want to keep in the final result.
       if(!edges1.empty())
         {
           try
@@ -8206,17 +8787,223 @@ void MEDCouplingUMesh::BuildIntersecting2DCellsFromEdges(double eps, const MEDCo
     }
 }
 
+void IKGeo2DInternalMapper2(INTERP_KERNEL::Node *n, const std::map<INTERP_KERNEL::Node *,int>& m, int forbVal0, int forbVal1, std::vector<int>& isect)
+{
+  std::map<INTERP_KERNEL::Node *,int>::const_iterator it(m.find(n));
+  if(it==m.end())
+    throw INTERP_KERNEL::Exception("Internal error in remapping !");
+  int v((*it).second);
+  if(v==forbVal0 || v==forbVal1)
+    return ;
+  if(std::find(isect.begin(),isect.end(),v)==isect.end())
+    isect.push_back(v);
+}
+
+bool IKGeo2DInternalMapper(const INTERP_KERNEL::ComposedEdge& c, const std::map<INTERP_KERNEL::Node *,int>& m, int forbVal0, int forbVal1, std::vector<int>& isect)
+{
+  int sz(c.size());
+  if(sz<=1)
+    return false;
+  bool presenceOfOn(false);
+  for(int i=0;i<sz;i++)
+    {
+      INTERP_KERNEL::ElementaryEdge *e(c[i]);
+      if(e->getLoc()!=INTERP_KERNEL::FULL_ON_1)
+        continue ;
+      IKGeo2DInternalMapper2(e->getStartNode(),m,forbVal0,forbVal1,isect);
+      IKGeo2DInternalMapper2(e->getEndNode(),m,forbVal0,forbVal1,isect);
+    }
+  return presenceOfOn;
+}
+
+/**
+ * This method split some of edges of 2D cells in \a this. The edges to be split are specified in \a subNodesInSeg and in \a subNodesInSegI using index storage mode.
+ * To do the work this method can optionnaly needs information about middle of subedges for quadratic cases if a minimal creation of new nodes is wanted.
+ * So this method try to reduce at most the number of new nodes. The only case that can lead this method to add nodes if a SEG3 is split without information of middle.
+ * \b WARNING : is returned value is different from 0 a call to MEDCouplingUMesh::mergeNodes is necessary to avoid to have a non conform mesh.
+ *
+ * \return int - the number of new nodes created (in most of cases 0).
+ * 
+ * \throw If \a this is not coherent.
+ * \throw If \a this has not spaceDim equal to 2.
+ * \throw If \a this has not meshDim equal to 2.
+ * \throw If some subcells needed to be split are orphan.
+ * \sa MEDCouplingUMesh::conformize2D
+ */
+int MEDCouplingUMesh::split2DCells(const DataArrayInt *desc, const DataArrayInt *descI, const DataArrayInt *subNodesInSeg, const DataArrayInt *subNodesInSegI, const DataArrayInt *midOpt, const DataArrayInt *midOptI)
+{
+  if(!desc || !descI || !subNodesInSeg || !subNodesInSegI)
+    throw INTERP_KERNEL::Exception("MEDCouplingUMesh::split2DCells : the 4 first arrays must be not null !");
+  desc->checkAllocated(); descI->checkAllocated(); subNodesInSeg->checkAllocated(); subNodesInSegI->checkAllocated();
+  if(getSpaceDimension()!=2 || getMeshDimension()!=2)
+    throw INTERP_KERNEL::Exception("MEDCouplingUMesh::split2DCells : This method only works for meshes with spaceDim=2 and meshDim=2 !");
+  if(midOpt==0 && midOptI==0)
+    {
+      split2DCellsLinear(desc,descI,subNodesInSeg,subNodesInSegI);
+      return 0;
+    }
+  else if(midOpt!=0 && midOptI!=0)
+    return split2DCellsQuadratic(desc,descI,subNodesInSeg,subNodesInSegI,midOpt,midOptI);
+  else
+    throw INTERP_KERNEL::Exception("MEDCouplingUMesh::split2DCells : middle parameters must be set to null for all or not null for all.");
+}
+
 /*!
- * This method is private and is the first step of Partition of 2D mesh (spaceDim==2 and meshDim==2).
+ * \b WARNING this method is \b potentially \b non \b const (if returned array is empty).
+ * \b WARNING this method lead to have a non geometric type sorted mesh (for MED file users) !
+ * This method performs a conformization of \b this. So if a edge in \a this can be split into entire edges in \a this this method
+ * will suppress such edges to use sub edges in \a this. So this method does not add nodes in \a this if merged edges have same nature each other (Linear,Quadratic).
+ * Whatever the returned value, this method does not alter the order of cells in \a this neither the orientation of cells.
+ * The modified cells if any are systematically declared as NORM_POLYGON or NORM_QPOLYG depending on the 
+ *
+ * This method expects that \b this has a meshDim equal 2 and spaceDim equal to 2 too.
+ * This method expects that all nodes in \a this are not closer than \a eps.
+ * If it is not the case you can invoke MEDCouplingUMesh::mergeNodes before calling this method.
  * 
+ * \param [in] eps the relative error to detect merged edges.
+ * \return DataArrayInt  * - The list of cellIds in \a this that have been subdivided. If empty, nothing changed in \a this (as if it were a const method). The array is a newly allocated array
+ *                           that the user is expected to deal with.
+ *
+ * \throw If \a this is not coherent.
+ * \throw If \a this has not spaceDim equal to 2.
+ * \throw If \a this has not meshDim equal to 2.
+ * \sa MEDCouplingUMesh::mergeNodes, MEDCouplingUMesh::split2DCells
+ */
+DataArrayInt *MEDCouplingUMesh::conformize2D(double eps)
+{
+  static const int SPACEDIM=2;
+  checkCoherency();
+  if(getSpaceDimension()!=2 || getMeshDimension()!=2)
+    throw INTERP_KERNEL::Exception("MEDCouplingUMesh::conformize2D : This method only works for meshes with spaceDim=2 and meshDim=2 !");
+  MEDCouplingAutoRefCountObjectPtr<DataArrayInt> desc1(DataArrayInt::New()),descIndx1(DataArrayInt::New()),revDesc1(DataArrayInt::New()),revDescIndx1(DataArrayInt::New());
+  MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> mDesc(buildDescendingConnectivity(desc1,descIndx1,revDesc1,revDescIndx1));
+  const int *c(mDesc->getNodalConnectivity()->getConstPointer()),*ci(mDesc->getNodalConnectivityIndex()->getConstPointer()),*rd(revDesc1->getConstPointer()),*rdi(revDescIndx1->getConstPointer());
+  MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> bboxArr(mDesc->getBoundingBoxForBBTree());
+  const double *bbox(bboxArr->begin()),*coords(getCoords()->begin());
+  int nCell(getNumberOfCells()),nDescCell(mDesc->getNumberOfCells());
+  std::vector< std::vector<int> > intersectEdge(nDescCell),overlapEdge(nDescCell);
+  std::vector<double> addCoo;
+  BBTree<SPACEDIM,int> myTree(bbox,0,0,nDescCell,-eps);
+  INTERP_KERNEL::QUADRATIC_PLANAR::_precision=eps;
+  INTERP_KERNEL::QUADRATIC_PLANAR::_arc_detection_precision=eps;
+  for(int i=0;i<nDescCell;i++)
+    {
+      std::vector<int> candidates;
+      myTree.getIntersectingElems(bbox+i*2*SPACEDIM,candidates);
+      for(std::vector<int>::const_iterator it=candidates.begin();it!=candidates.end();it++)
+        if(*it>i)
+          {
+            std::map<INTERP_KERNEL::Node *,int> m;
+            INTERP_KERNEL::Edge *e1(MEDCouplingUMeshBuildQPFromEdge2((INTERP_KERNEL::NormalizedCellType)c[ci[i]],c+ci[i]+1,coords,m)),
+              *e2(MEDCouplingUMeshBuildQPFromEdge2((INTERP_KERNEL::NormalizedCellType)c[ci[*it]],c+ci[*it]+1,coords,m));
+            INTERP_KERNEL::MergePoints merge;
+            INTERP_KERNEL::QuadraticPolygon c1,c2;
+            e1->intersectWith(e2,merge,c1,c2);
+            e1->decrRef(); e2->decrRef();
+            if(IKGeo2DInternalMapper(c1,m,c[ci[i]+1],c[ci[i]+2],intersectEdge[i]))
+              overlapEdge[i].push_back(*it);
+            if(IKGeo2DInternalMapper(c2,m,c[ci[*it]+1],c[ci[*it]+2],intersectEdge[*it]))
+              overlapEdge[*it].push_back(i);
+            for(std::map<INTERP_KERNEL::Node *,int>::const_iterator it2=m.begin();it2!=m.end();it2++)
+              (*it2).first->decrRef();
+          }
+    }
+  // splitting done. sort intersect point in intersectEdge.
+  std::vector< std::vector<int> > middle(nDescCell);
+  int nbOf2DCellsToBeSplit(0);
+  bool middleNeedsToBeUsed(false);
+  std::vector<bool> cells2DToTreat(nDescCell,false);
+  for(int i=0;i<nDescCell;i++)
+    {
+      std::vector<int>& isect(intersectEdge[i]);
+      int sz((int)isect.size());
+      if(sz>1)
+        {
+          std::map<INTERP_KERNEL::Node *,int> m;
+          INTERP_KERNEL::Edge *e(MEDCouplingUMeshBuildQPFromEdge2((INTERP_KERNEL::NormalizedCellType)c[ci[i]],c+ci[i]+1,coords,m));
+          e->sortSubNodesAbs(coords,isect);
+          e->decrRef();
+          for(std::map<INTERP_KERNEL::Node *,int>::const_iterator it2=m.begin();it2!=m.end();it2++)
+            (*it2).first->decrRef();
+        }
+      if(sz!=0)
+        {
+          int idx0(rdi[i]),idx1(rdi[i+1]);
+          if(idx1-idx0!=1)
+            throw INTERP_KERNEL::Exception("MEDCouplingUMesh::conformize2D : internal error #0 !");
+          if(!cells2DToTreat[rd[idx0]])
+            {
+              cells2DToTreat[rd[idx0]]=true;
+              nbOf2DCellsToBeSplit++;
+            }
+          // try to reuse at most eventual 'middle' of SEG3
+          std::vector<int>& mid(middle[i]);
+          mid.resize(sz+1,-1);
+          if((INTERP_KERNEL::NormalizedCellType)c[ci[i]]==INTERP_KERNEL::NORM_SEG3)
+            {
+              middleNeedsToBeUsed=true;
+              const std::vector<int>& candidates(overlapEdge[i]);
+              std::vector<int> trueCandidates;
+              for(std::vector<int>::const_iterator itc=candidates.begin();itc!=candidates.end();itc++)
+                if((INTERP_KERNEL::NormalizedCellType)c[ci[*itc]]==INTERP_KERNEL::NORM_SEG3)
+                  trueCandidates.push_back(*itc);
+              int stNode(c[ci[i]+1]),endNode(isect[0]);
+              for(int j=0;j<sz+1;j++)
+                {
+                  for(std::vector<int>::const_iterator itc=trueCandidates.begin();itc!=trueCandidates.end();itc++)
+                    {
+                      int tmpSt(c[ci[*itc]+1]),tmpEnd(c[ci[*itc]+2]);
+                      if((tmpSt==stNode && tmpEnd==endNode) || (tmpSt==endNode && tmpEnd==stNode))
+                        { mid[j]=*itc; break; }
+                    }
+                  stNode=endNode;
+                  endNode=j<sz-1?isect[j+1]:c[ci[i]+2];
+                }
+            }
+        }
+    }
+  MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret(DataArrayInt::New()),notRet(DataArrayInt::New()); ret->alloc(nbOf2DCellsToBeSplit,1);
+  if(nbOf2DCellsToBeSplit==0)
+    return ret.retn();
+  //
+  int *retPtr(ret->getPointer());
+  for(int i=0;i<nCell;i++)
+    if(cells2DToTreat[i])
+      *retPtr++=i;
+  //
+  MEDCouplingAutoRefCountObjectPtr<DataArrayInt> mSafe,nSafe,oSafe,pSafe,qSafe,rSafe;
+  DataArrayInt *m(0),*n(0),*o(0),*p(0),*q(0),*r(0);
+  MEDCouplingUMesh::ExtractFromIndexedArrays(ret->begin(),ret->end(),desc1,descIndx1,m,n); mSafe=m; nSafe=n;
+  DataArrayInt::PutIntoToSkylineFrmt(intersectEdge,o,p); oSafe=o; pSafe=p;
+  if(middleNeedsToBeUsed)
+    { DataArrayInt::PutIntoToSkylineFrmt(middle,q,r); qSafe=q; rSafe=r; }
+  MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> modif(static_cast<MEDCouplingUMesh *>(buildPartOfMySelf(ret->begin(),ret->end(),true)));
+  int nbOfNodesCreated(modif->split2DCells(mSafe,nSafe,oSafe,pSafe,qSafe,rSafe));
+  setCoords(modif->getCoords());//if nbOfNodesCreated==0 modif and this have the same coordinates pointer so this line has no effect. But for quadratic cases this line is important.
+  setPartOfMySelf(ret->begin(),ret->end(),*modif);
+  {
+    bool areNodesMerged; int newNbOfNodes;
+    if(nbOfNodesCreated!=0)
+      MEDCouplingAutoRefCountObjectPtr<DataArrayInt> tmp(mergeNodes(eps,areNodesMerged,newNbOfNodes));
+  }
+  return ret.retn();
+}
+
+/*!
+ * This method is private and is the first step of Partition of 2D mesh (spaceDim==2 and meshDim==2).
+ * It builds the descending connectivity of the two meshes, and then using a binary tree
+ * it computes the edge intersections. This results in new points being created : they're stored in addCoo.
+ * Documentation about parameters  colinear2 and subDiv2 can be found in method QuadraticPolygon::splitAbs().
  */
 void MEDCouplingUMesh::IntersectDescending2DMeshes(const MEDCouplingUMesh *m1, const MEDCouplingUMesh *m2, double eps,
                                                    std::vector< std::vector<int> >& intersectEdge1, std::vector< std::vector<int> >& colinear2, std::vector< std::vector<int> >& subDiv2,
                                                    MEDCouplingUMesh *& m1Desc, DataArrayInt *&desc1, DataArrayInt *&descIndx1, DataArrayInt *&revDesc1, DataArrayInt *&revDescIndx1,
-                                                   MEDCouplingUMesh *& m2Desc, DataArrayInt *&desc2, DataArrayInt *&descIndx2, DataArrayInt *&revDesc2, DataArrayInt *&revDescIndx2,
-                                                   std::vector<double>& addCoo) throw(INTERP_KERNEL::Exception)
+                                                   std::vector<double>& addCoo,
+                                                   MEDCouplingUMesh *& m2Desc, DataArrayInt *&desc2, DataArrayInt *&descIndx2, DataArrayInt *&revDesc2, DataArrayInt *&revDescIndx2)
+                                                   throw(INTERP_KERNEL::Exception)
 {
   static const int SPACEDIM=2;
+  // Build desc connectivity
   desc1=DataArrayInt::New(); descIndx1=DataArrayInt::New(); revDesc1=DataArrayInt::New(); revDescIndx1=DataArrayInt::New();
   desc2=DataArrayInt::New();
   descIndx2=DataArrayInt::New();
@@ -8229,28 +9016,42 @@ 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);
-  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);
+
+  // Build BB tree of all edges in the tool mesh (second mesh)
+  MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> bbox1Arr(m1Desc->getBoundingBoxForBBTree()),bbox2Arr(m2Desc->getBoundingBoxForBBTree());
+  const double *bbox1(bbox1Arr->begin()),*bbox2(bbox2Arr->begin());
+  int nDescCell1=m1Desc->getNumberOfCells();
+  int nDescCell2=m2Desc->getNumberOfCells();
+  intersectEdge1.resize(nDescCell1);
+  colinear2.resize(nDescCell2);
+  subDiv2.resize(nDescCell2);
+  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++)
+  for(int i=0;i<nDescCell1;i++)  // for all edges in the first mesh
     {
-      std::vector<int> candidates2;
-      myTree.getIntersectingElems(&bbox1[i*2*SPACEDIM],candidates2);
-      if(!candidates2.empty())
+      std::vector<int> candidates2; // edges of mesh2 candidate for intersection
+      myTree.getIntersectingElems(bbox1+i*2*SPACEDIM,candidates2);
+      if(!candidates2.empty()) // candidates2 holds edges from the second mesh potentially intersecting current edge i in mesh1
         {
           std::map<INTERP_KERNEL::Node *,int> map1,map2;
+          // pol2 is not necessarily a closed polygon: just a set of (quadratic) edges (same as candidates2) in the Geometric DS format
           INTERP_KERNEL::QuadraticPolygon *pol2=MEDCouplingUMeshBuildQPFromMesh(m2Desc,candidates2,map2);
           candidates1[0]=i;
           INTERP_KERNEL::QuadraticPolygon *pol1=MEDCouplingUMeshBuildQPFromMesh(m1Desc,candidates1,map1);
+          // This following part is to avoid that some removed nodes (for example due to a merge between pol1 and pol2) are replaced by a newly created one
+          // This trick guarantees that Node * are discriminant (i.e. form a unique identifier)
+          std::set<INTERP_KERNEL::Node *> nodes;
+          pol1->getAllNodes(nodes); pol2->getAllNodes(nodes);
+          std::size_t szz(nodes.size());
+          std::vector< MEDCouplingAutoRefCountObjectPtr<INTERP_KERNEL::Node> > nodesSafe(szz);
+          std::set<INTERP_KERNEL::Node *>::const_iterator itt(nodes.begin());
+          for(std::size_t iii=0;iii<szz;iii++,itt++)
+            { (*itt)->incrRef(); nodesSafe[iii]=*itt; }
+          // end of protection
+          // Performs egde cutting:
           pol1->splitAbs(*pol2,map1,map2,offset1,offset2,candidates2,intersectEdge1[i],i,colinear2,subDiv2,addCoo);
           delete pol2;
           delete pol1;
@@ -8267,13 +9068,18 @@ void MEDCouplingUMesh::IntersectDescending2DMeshes(const MEDCouplingUMesh *m1, c
  * This method has 4 inputs :
  *  - a mesh 'm1' with meshDim==1 and a SpaceDim==2
  *  - a mesh 'm2' with meshDim==1 and a SpaceDim==2
- *  - subDiv of size 'm2->getNumberOfCells()' that lists for each seg cell in 'm' the splitting node ids in randomly sorted.
- * The aim of this method is to sort the splitting nodes, if any, and to put in 'intersectEdge' output paramter based on edges of mesh 'm2'
- * \param m1 is expected to be a mesh of meshDimension equal to 1 and spaceDim equal to 2. No check of that is performed by this method. Only present for its coords in case of 'subDiv' shares some nodes of 'm1'
+ *  - subDiv of size 'm2->getNumberOfCells()' that lists for each seg cell in 'm' the splitting node ids randomly sorted.
+ * The aim of this method is to sort the splitting nodes, if any, and to put them in 'intersectEdge' output parameter based on edges of mesh 'm2'
+ * Nodes end up lying consecutively on a cutted edge.
+ * \param m1 is expected to be a mesh of meshDimension equal to 1 and spaceDim equal to 2. No check of that is performed by this method.
+ * (Only present for its coords in case of 'subDiv' shares some nodes of 'm1')
  * \param m2 is expected to be a mesh of meshDimension equal to 1 and spaceDim equal to 2. No check of that is performed by this method.
- * \param addCoo input parameter with additionnal nodes linked to intersection of the 2 meshes.
+ * \param addCoo input parameter with additional nodes linked to intersection of the 2 meshes.
+ * \param[out] intersectEdge the same content as subDiv, but correclty oriented.
  */
-void MEDCouplingUMesh::BuildIntersectEdges(const MEDCouplingUMesh *m1, const MEDCouplingUMesh *m2, const std::vector<double>& addCoo, const std::vector< std::vector<int> >& subDiv, std::vector< std::vector<int> >& intersectEdge) throw(INTERP_KERNEL::Exception)
+void MEDCouplingUMesh::BuildIntersectEdges(const MEDCouplingUMesh *m1, const MEDCouplingUMesh *m2,
+      const std::vector<double>& addCoo,
+      const std::vector< std::vector<int> >& subDiv, std::vector< std::vector<int> >& intersectEdge)
 {
   int offset1=m1->getNumberOfNodes();
   int ncell=m2->getNumberOfCells();
@@ -8490,7 +9296,7 @@ void MEDCouplingUMesh::assemblyForSplitFrom3DSurf(const std::vector< std::pair<i
  * 
  * \return false if the input connectivity represents already the convex hull, true if the input cell needs to be reordered.
  */
-bool MEDCouplingUMesh::BuildConvexEnvelopOf2DCellJarvis(const double *coords, const int *nodalConnBg, const int *nodalConnEnd, DataArrayInt *nodalConnecOut) throw(INTERP_KERNEL::Exception)
+bool MEDCouplingUMesh::BuildConvexEnvelopOf2DCellJarvis(const double *coords, const int *nodalConnBg, const int *nodalConnEnd, DataArrayInt *nodalConnecOut)
 {
   std::size_t sz=std::distance(nodalConnBg,nodalConnEnd);
   if(sz>=4)
@@ -8589,7 +9395,7 @@ bool MEDCouplingUMesh::BuildConvexEnvelopOf2DCellJarvis(const double *coords, co
  * \param [in] offsetForRemoval (by default 0) offset so that for each i in [0,arrIndx->getNumberOfTuples()-1) removal process will be performed in the following range [arr+arrIndx[i]+offsetForRemoval,arr+arr[i+1])
  * \return true if \b arr and \b arrIndx have been modified, false if not.
  */
-bool MEDCouplingUMesh::RemoveIdsFromIndexedArrays(const int *idsToRemoveBg, const int *idsToRemoveEnd, DataArrayInt *arr, DataArrayInt *arrIndx, int offsetForRemoval) throw(INTERP_KERNEL::Exception)
+bool MEDCouplingUMesh::RemoveIdsFromIndexedArrays(const int *idsToRemoveBg, const int *idsToRemoveEnd, DataArrayInt *arr, DataArrayInt *arrIndx, int offsetForRemoval)
 {
   if(!arrIndx || !arr)
     throw INTERP_KERNEL::Exception("MEDCouplingUMesh::RemoveIdsFromIndexedArrays : some input arrays are empty !");
@@ -8730,7 +9536,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];
@@ -8751,7 +9557,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);
@@ -8890,7 +9696,7 @@ void MEDCouplingUMesh::SetPartOfIndexedArraysSameIdx(const int *idsOfSelectBg, c
  * \return a newly allocated DataArray that stores all ids fetched by the gradually spread process.
  * \sa MEDCouplingUMesh::ComputeSpreadZoneGraduallyFromSeed, MEDCouplingUMesh::partitionBySpreadZone
  */
-DataArrayInt *MEDCouplingUMesh::ComputeSpreadZoneGradually(const DataArrayInt *arrIn, const DataArrayInt *arrIndxIn) throw(INTERP_KERNEL::Exception)
+DataArrayInt *MEDCouplingUMesh::ComputeSpreadZoneGradually(const DataArrayInt *arrIn, const DataArrayInt *arrIndxIn)
 {
   int seed=0,nbOfDepthPeelingPerformed=0;
   return ComputeSpreadZoneGraduallyFromSeed(&seed,&seed+1,arrIn,arrIndxIn,-1,nbOfDepthPeelingPerformed);
@@ -8912,7 +9718,7 @@ DataArrayInt *MEDCouplingUMesh::ComputeSpreadZoneGradually(const DataArrayInt *a
  * \return a newly allocated DataArray that stores all ids fetched by the gradually spread process.
  * \sa MEDCouplingUMesh::partitionBySpreadZone
  */
-DataArrayInt *MEDCouplingUMesh::ComputeSpreadZoneGraduallyFromSeed(const int *seedBg, const int *seedEnd, const DataArrayInt *arrIn, const DataArrayInt *arrIndxIn, int nbOfDepthPeeling, int& nbOfDepthPeelingPerformed) throw(INTERP_KERNEL::Exception)
+DataArrayInt *MEDCouplingUMesh::ComputeSpreadZoneGraduallyFromSeed(const int *seedBg, const int *seedEnd, const DataArrayInt *arrIn, const DataArrayInt *arrIndxIn, int nbOfDepthPeeling, int& nbOfDepthPeelingPerformed)
 {
   nbOfDepthPeelingPerformed=0;
   if(!arrIndxIn)
@@ -8928,7 +9734,7 @@ DataArrayInt *MEDCouplingUMesh::ComputeSpreadZoneGraduallyFromSeed(const int *se
   return ComputeSpreadZoneGraduallyFromSeedAlg(fetched,seedBg,seedEnd,arrIn,arrIndxIn,nbOfDepthPeeling,nbOfDepthPeelingPerformed);
 }
 
-DataArrayInt *MEDCouplingUMesh::ComputeSpreadZoneGraduallyFromSeedAlg(std::vector<bool>& fetched, const int *seedBg, const int *seedEnd, const DataArrayInt *arrIn, const DataArrayInt *arrIndxIn, int nbOfDepthPeeling, int& nbOfDepthPeelingPerformed) throw(INTERP_KERNEL::Exception)
+DataArrayInt *MEDCouplingUMesh::ComputeSpreadZoneGraduallyFromSeedAlg(std::vector<bool>& fetched, const int *seedBg, const int *seedEnd, const DataArrayInt *arrIn, const DataArrayInt *arrIndxIn, int nbOfDepthPeeling, int& nbOfDepthPeelingPerformed)
 {
   nbOfDepthPeelingPerformed=0;
   if(!seedBg || !seedEnd || !arrIn || !arrIndxIn)
@@ -9092,7 +9898,7 @@ void MEDCouplingUMesh::SetPartOfIndexedArraysSameIdx2(int start, int end, int st
  * 
  * \return a newly allocated mesh lying on the same coords than \b this with same meshdimension than \b this.
  */
-MEDCouplingUMesh *MEDCouplingUMesh::buildSpreadZonesWithPoly() const throw(INTERP_KERNEL::Exception)
+MEDCouplingUMesh *MEDCouplingUMesh::buildSpreadZonesWithPoly() const
 {
   checkFullyDefined();
   int mdim=getMeshDimension();
@@ -9134,9 +9940,8 @@ MEDCouplingUMesh *MEDCouplingUMesh::buildSpreadZonesWithPoly() const throw(INTER
  * This method only needs a well defined connectivity. Coordinates are not considered here.
  * This method returns a vector of \b newly allocated arrays that the caller has to deal with.
  */
-std::vector<DataArrayInt *> MEDCouplingUMesh::partitionBySpreadZone() const throw(INTERP_KERNEL::Exception)
+std::vector<DataArrayInt *> MEDCouplingUMesh::partitionBySpreadZone() const
 {
-  //#if 0
   int nbOfCellsCur=getNumberOfCells();
   std::vector<DataArrayInt *> ret;
   if(nbOfCellsCur<=0)
@@ -9156,36 +9961,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
 }
 
 /*!
@@ -9196,7 +9971,7 @@ std::vector<DataArrayInt *> MEDCouplingUMesh::partitionBySpreadZone() const thro
  * \return a newly allocated DataArrayInt to be managed by the caller.
  * \throw In case of \a code has not the right format (typically of size 3*n)
  */
-DataArrayInt *MEDCouplingUMesh::ComputeRangesFromTypeDistribution(const std::vector<int>& code) throw(INTERP_KERNEL::Exception)
+DataArrayInt *MEDCouplingUMesh::ComputeRangesFromTypeDistribution(const std::vector<int>& code)
 {
   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret=DataArrayInt::New();
   std::size_t nb=code.size()/3;
@@ -9212,6 +9987,191 @@ 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.
+ *
+ * \warning This method operates on each cells in this independantly ! So it can leads to non conform mesh in returned value ! If you expect to have a conform mesh in output
+ * the policy PLANAR_FACE_6 should be used on a mesh sorted with MEDCoupling1SGTUMesh::sortHexa8EachOther.
+ * 
+ * \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::sortHexa8EachOther
+ */
+MEDCoupling1SGTUMesh *MEDCouplingUMesh::tetrahedrize(int policy, DataArrayInt *& n2oCells, int& nbOfAdditionalPoints) const
+{
+  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(),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();
+}
+
+/*!
+ * It is the linear part of MEDCouplingUMesh::split2DCells. Here no additionnal nodes will be added in \b this. So coordinates pointer remain unchanged (is not even touch). 
+ *
+ * \sa MEDCouplingUMesh::split2DCells
+ */
+void MEDCouplingUMesh::split2DCellsLinear(const DataArrayInt *desc, const DataArrayInt *descI, const DataArrayInt *subNodesInSeg, const DataArrayInt *subNodesInSegI)
+{
+  checkConnectivityFullyDefined();
+  int ncells(getNumberOfCells()),lgthToReach(getMeshLength()+subNodesInSeg->getNumberOfTuples());
+  MEDCouplingAutoRefCountObjectPtr<DataArrayInt> c(DataArrayInt::New()); c->alloc((std::size_t)lgthToReach);
+  const int *subPtr(subNodesInSeg->begin()),*subIPtr(subNodesInSegI->begin()),*descPtr(desc->begin()),*descIPtr(descI->begin()),*oldConn(getNodalConnectivity()->begin());
+  int *cPtr(c->getPointer()),*ciPtr(getNodalConnectivityIndex()->getPointer());
+  int prevPosOfCi(ciPtr[0]);
+  for(int i=0;i<ncells;i++,ciPtr++,descIPtr++)
+    {
+      int offset(descIPtr[0]),sz(descIPtr[1]-descIPtr[0]),deltaSz(0);
+      *cPtr++=(int)INTERP_KERNEL::NORM_POLYGON; *cPtr++=oldConn[prevPosOfCi+1];
+      for(int j=0;j<sz;j++)
+        {
+          int offset2(subIPtr[descPtr[offset+j]]),sz2(subIPtr[descPtr[offset+j]+1]-subIPtr[descPtr[offset+j]]);
+          for(int k=0;k<sz2;k++)
+            *cPtr++=subPtr[offset2+k];
+          if(j!=sz-1)
+            *cPtr++=oldConn[prevPosOfCi+j+2];
+          deltaSz+=sz2;
+        }
+      prevPosOfCi=ciPtr[1];
+      ciPtr[1]=ciPtr[0]+1+sz+deltaSz;//sz==old nb of nodes because (nb of subedges=nb of nodes for polygons)
+    }
+  if(c->end()!=cPtr)
+    throw INTERP_KERNEL::Exception("MEDCouplingUMesh::split2DCellsLinear : Some of edges to be split are orphan !");
+  _nodal_connec->decrRef();
+  _nodal_connec=c.retn(); _types.clear(); _types.insert(INTERP_KERNEL::NORM_POLYGON);
+}
+
+int internalAddPoint(const INTERP_KERNEL::Edge *e, int id, const double *coo, int startId, int endId, DataArrayDouble& addCoo, int& nodesCnter)
+{
+  if(id!=-1)
+    return id;
+  else
+    {
+      int ret(nodesCnter++);
+      double newPt[2];
+      e->getMiddleOfPoints(coo+2*startId,coo+2*endId,newPt);
+      addCoo.insertAtTheEnd(newPt,newPt+2);
+      return ret;
+    }
+}
+
+/*!
+ * It is the quadratic part of MEDCouplingUMesh::split2DCells. Here some additionnal nodes can be added at the end of coordinates array object.
+ *
+ * \return  int - the number of new nodes created.
+ * \sa MEDCouplingUMesh::split2DCells
+ */
+int MEDCouplingUMesh::split2DCellsQuadratic(const DataArrayInt *desc, const DataArrayInt *descI, const DataArrayInt *subNodesInSeg, const DataArrayInt *subNodesInSegI, const DataArrayInt *mid, const DataArrayInt *midI)
+{
+  checkCoherency();
+  int ncells(getNumberOfCells()),lgthToReach(getMeshLength()+2*subNodesInSeg->getNumberOfTuples()),nodesCnt(getNumberOfNodes());
+  MEDCouplingAutoRefCountObjectPtr<DataArrayInt> c(DataArrayInt::New()); c->alloc((std::size_t)lgthToReach);
+  MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> addCoo(DataArrayDouble::New()); addCoo->alloc(0,1);
+  const int *subPtr(subNodesInSeg->begin()),*subIPtr(subNodesInSegI->begin()),*descPtr(desc->begin()),*descIPtr(descI->begin()),*oldConn(getNodalConnectivity()->begin());
+  const int *midPtr(mid->begin()),*midIPtr(midI->begin());
+  const double *oldCoordsPtr(getCoords()->begin());
+  int *cPtr(c->getPointer()),*ciPtr(getNodalConnectivityIndex()->getPointer());
+  int prevPosOfCi(ciPtr[0]);
+  for(int i=0;i<ncells;i++,ciPtr++,descIPtr++)
+    {
+      int offset(descIPtr[0]),sz(descIPtr[1]-descIPtr[0]),deltaSz(sz);
+      for(int j=0;j<sz;j++)
+        { int sz2(subIPtr[descPtr[offset+j]+1]-subIPtr[descPtr[offset+j]]); deltaSz+=sz2; }
+      *cPtr++=(int)INTERP_KERNEL::NORM_QPOLYG; cPtr[0]=oldConn[prevPosOfCi+1];
+      for(int j=0;j<sz;j++)//loop over subedges of oldConn
+        {
+          int offset2(subIPtr[descPtr[offset+j]]),sz2(subIPtr[descPtr[offset+j]+1]-subIPtr[descPtr[offset+j]]),offset3(midIPtr[descPtr[offset+j]]);
+          if(sz2==0)
+            {
+              if(j<sz-1)
+                cPtr[1]=oldConn[prevPosOfCi+2+j];
+              cPtr[deltaSz]=oldConn[prevPosOfCi+1+j+sz]; cPtr++;
+              continue;
+            }
+          std::vector<INTERP_KERNEL::Node *> ns(3);
+          ns[0]=new INTERP_KERNEL::Node(oldCoordsPtr[2*oldConn[prevPosOfCi+1+j]],oldCoordsPtr[2*oldConn[prevPosOfCi+1+j]+1]);
+          ns[1]=new INTERP_KERNEL::Node(oldCoordsPtr[2*oldConn[prevPosOfCi+1+(1+j)%sz]],oldCoordsPtr[2*oldConn[prevPosOfCi+1+(1+j)%sz]+1]);
+          ns[2]=new INTERP_KERNEL::Node(oldCoordsPtr[2*oldConn[prevPosOfCi+1+sz+j]],oldCoordsPtr[2*oldConn[prevPosOfCi+1+sz+j]+1]);
+          MEDCouplingAutoRefCountObjectPtr<INTERP_KERNEL::Edge> e(INTERP_KERNEL::QuadraticPolygon::BuildArcCircleEdge(ns));
+          for(int k=0;k<sz2;k++)//loop over subsplit of current subedge
+            {
+              cPtr[1]=subPtr[offset2+k];
+              cPtr[deltaSz]=internalAddPoint(e,midPtr[offset3+k],oldCoordsPtr,cPtr[0],cPtr[1],*addCoo,nodesCnt); cPtr++;
+            }
+          int tmpEnd(oldConn[prevPosOfCi+1+(j+1)%sz]);
+          if(j!=sz-1)
+            { cPtr[1]=tmpEnd; }
+          cPtr[deltaSz]=internalAddPoint(e,midPtr[offset3+sz2],oldCoordsPtr,cPtr[0],tmpEnd,*addCoo,nodesCnt); cPtr++;
+        }
+      prevPosOfCi=ciPtr[1]; cPtr+=deltaSz;
+      ciPtr[1]=ciPtr[0]+1+2*deltaSz;//sz==old nb of nodes because (nb of subedges=nb of nodes for polygons)
+    }
+  if(c->end()!=cPtr)
+    throw INTERP_KERNEL::Exception("MEDCouplingUMesh::split2DCellsQuadratic : Some of edges to be split are orphan !");
+  _nodal_connec->decrRef();
+  _nodal_connec=c.retn(); _types.clear(); _types.insert(INTERP_KERNEL::NORM_QPOLYG);
+  addCoo->rearrange(2);
+  MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> coo(DataArrayDouble::Aggregate(getCoords(),addCoo));//info are copied from getCoords() by using Aggregate
+  setCoords(coo);
+  return addCoo->getNumberOfTuples();
+}
+
 MEDCouplingUMeshCellIterator::MEDCouplingUMeshCellIterator(MEDCouplingUMesh *mesh):_mesh(mesh),_cell(new MEDCouplingUMeshCell(mesh)),
                                                                                    _own_cell(true),_cell_id(-1),_nb_cell(0)
 {