X-Git-Url: http://git.salome-platform.org/gitweb/?a=blobdiff_plain;f=src%2FMEDCoupling%2FMEDCouplingUMesh.cxx;h=b19c250b7601d6f513a291f02bfec95cc6ff463c;hb=ac1df6b0ba8b337555fb39610c89f678d889580d;hp=33f9e99930669be941db2d91b37d2dcca172dea0;hpb=51fff402deffba5b8c47ae68f4c26f1924b700bb;p=tools%2Fmedcoupling.git

diff --git a/src/MEDCoupling/MEDCouplingUMesh.cxx b/src/MEDCoupling/MEDCouplingUMesh.cxx
index 33f9e9993..b19c250b7 100644
--- a/src/MEDCoupling/MEDCouplingUMesh.cxx
+++ b/src/MEDCoupling/MEDCouplingUMesh.cxx
@@ -16,7 +16,7 @@
 //
 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
 //
-// Author : Anthony Geay (CEA/DEN)
+// Author : Anthony Geay (EDF R&D)
 
 #include "MEDCouplingUMesh.hxx"
 #include "MEDCouplingCMesh.hxx"
@@ -56,8 +56,8 @@ using namespace MEDCoupling;
 double MEDCouplingUMesh::EPS_FOR_POLYH_ORIENTATION=1.e-14;
 
 /// @cond INTERNAL
-const INTERP_KERNEL::NormalizedCellType MEDCouplingUMesh::MEDMEM_ORDER[N_MEDMEM_ORDER] = { INTERP_KERNEL::NORM_POINT1, INTERP_KERNEL::NORM_SEG2, INTERP_KERNEL::NORM_SEG3, INTERP_KERNEL::NORM_SEG4, INTERP_KERNEL::NORM_POLYL, INTERP_KERNEL::NORM_TRI3, INTERP_KERNEL::NORM_QUAD4, INTERP_KERNEL::NORM_TRI6, INTERP_KERNEL::NORM_TRI7, INTERP_KERNEL::NORM_QUAD8, INTERP_KERNEL::NORM_QUAD9, INTERP_KERNEL::NORM_POLYGON, INTERP_KERNEL::NORM_QPOLYG, INTERP_KERNEL::NORM_TETRA4, INTERP_KERNEL::NORM_PYRA5, INTERP_KERNEL::NORM_PENTA6, INTERP_KERNEL::NORM_HEXA8, INTERP_KERNEL::NORM_HEXGP12, INTERP_KERNEL::NORM_TETRA10, INTERP_KERNEL::NORM_PYRA13, INTERP_KERNEL::NORM_PENTA15, INTERP_KERNEL::NORM_HEXA20, INTERP_KERNEL::NORM_HEXA27, INTERP_KERNEL::NORM_POLYHED };
-const int MEDCouplingUMesh::MEDCOUPLING2VTKTYPETRADUCER[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};
+const INTERP_KERNEL::NormalizedCellType MEDCouplingUMesh::MEDMEM_ORDER[N_MEDMEM_ORDER] = { INTERP_KERNEL::NORM_POINT1, INTERP_KERNEL::NORM_SEG2, INTERP_KERNEL::NORM_SEG3, INTERP_KERNEL::NORM_SEG4, INTERP_KERNEL::NORM_POLYL, INTERP_KERNEL::NORM_TRI3, INTERP_KERNEL::NORM_QUAD4, INTERP_KERNEL::NORM_TRI6, INTERP_KERNEL::NORM_TRI7, INTERP_KERNEL::NORM_QUAD8, INTERP_KERNEL::NORM_QUAD9, INTERP_KERNEL::NORM_POLYGON, INTERP_KERNEL::NORM_QPOLYG, INTERP_KERNEL::NORM_TETRA4, INTERP_KERNEL::NORM_PYRA5, INTERP_KERNEL::NORM_PENTA6, INTERP_KERNEL::NORM_HEXA8, INTERP_KERNEL::NORM_HEXGP12, INTERP_KERNEL::NORM_TETRA10, INTERP_KERNEL::NORM_PYRA13, INTERP_KERNEL::NORM_PENTA15, INTERP_KERNEL::NORM_PENTA18, INTERP_KERNEL::NORM_HEXA20, INTERP_KERNEL::NORM_HEXA27, INTERP_KERNEL::NORM_POLYHED };
+const int MEDCouplingUMesh::MEDCOUPLING2VTKTYPETRADUCER[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,32,-1,25,42,36,4};
 /// @endcond
 
 MEDCouplingUMesh *MEDCouplingUMesh::New()
@@ -301,9 +301,10 @@ void MEDCouplingUMesh::setMeshDimension(int meshDim)
 }
 
 /*!
- * Allocates memory to store an estimation of the given number of cells. The closer is the estimation to the number of cells effectively inserted,
- * the less will the library need to reallocate memory. If the number of cells to be inserted is not known simply put 0 to this parameter.
- * If a nodal connectivity previouly existed before the call of this method, it will be reset.
+ * Allocates memory to store an estimation of the given number of cells. 
+ * The closer the estimation to the number of cells effectively inserted, the less need the library requires
+ * to reallocate memory. If the number of cells to be inserted is not known simply assign 0 to this parameter.
+ * If a nodal connectivity previously existed before the call of this method, it will be reset.
  *
  *  \param [in] nbOfCells - estimation of the number of cell \a this mesh will contain.
  *
@@ -975,7 +976,12 @@ void MEDCouplingUMesh::computeEnlargedNeighborsOfNodes(MCAuto<DataArrayInt> &nei
   {
     int *neighIdx(neighborsIdx->getPointer());
     for(std::vector< std::set<int> >::const_iterator it=st0.begin();it!=st0.end();it++,neighIdx++)
-      neighIdx[1]=neighIdx[0]+(*it).size()-1;
+      {
+        if ((*it).empty())
+          neighIdx[1]=neighIdx[0];
+        else
+          neighIdx[1]=neighIdx[0]+(*it).size()-1;
+      }
   }
   neighbors=DataArrayInt::New(); neighbors->alloc(neighborsIdx->back(),1);
   {
@@ -1308,12 +1314,14 @@ void MEDCouplingUMesh::simplifyPolyhedra(double eps)
   connINew->alloc(nbOfCells+1,1);
   int *connINewPtr=connINew->getPointer(); *connINewPtr++=0;
   MCAuto<DataArrayInt> connNew=DataArrayInt::New(); connNew->alloc(0,1);
+  MCAuto<DataArrayInt> E_Fi(DataArrayInt::New()), E_F(DataArrayInt::New()), F_Ei(DataArrayInt::New()), F_E(DataArrayInt::New());
+  MCAuto<MEDCouplingUMesh> m_faces(buildDescendingConnectivity(E_F, E_Fi, F_E, F_Ei));
   bool changed=false;
   for(int i=0;i<nbOfCells;i++,connINewPtr++)
     {
       if(conn[index[i]]==(int)INTERP_KERNEL::NORM_POLYHED)
         {
-          SimplifyPolyhedronCell(eps,coords,conn+index[i],conn+index[i+1],connNew);
+          SimplifyPolyhedronCell(eps,coords, i,connNew, m_faces, E_Fi, E_F, F_Ei, F_E);
           changed=true;
         }
       else
@@ -1895,7 +1903,7 @@ MEDCouplingUMesh *MEDCouplingUMesh::mergeMyselfWithOnSameCoords(const MEDCouplin
 /*!
  * Build a sub part of \b this lying or not on the same coordinates than \b this (regarding value of \b keepCoords).
  * By default coordinates are kept. This method is close to MEDCouplingUMesh::buildPartOfMySelf except that here input
- * cellIds is not given explicitely but by a range python like.
+ * cellIds is not given explicitly but by a range python like.
  * 
  * \param start
  * \param end
@@ -1984,19 +1992,19 @@ void MEDCouplingUMesh::setPartOfMySelf(const int *cellIdsBg, const int *cellIdsE
       oss << ", whereas other mesh dimension is set equal to " << otherOnSameCoordsThanThis.getMeshDimension() << " !";
       throw INTERP_KERNEL::Exception(oss.str());
     }
-  int nbOfCellsToModify=(int)std::distance(cellIdsBg,cellIdsEnd);
+  std::size_t nbOfCellsToModify(std::distance(cellIdsBg,cellIdsEnd));
   if(nbOfCellsToModify!=otherOnSameCoordsThanThis.getNumberOfCells())
     {
       std::ostringstream oss; oss << "MEDCouplingUMesh::setPartOfMySelf : cells ids length (" <<  nbOfCellsToModify << ") do not match the number of cells of other mesh (" << otherOnSameCoordsThanThis.getNumberOfCells() << ") !";
       throw INTERP_KERNEL::Exception(oss.str());
     }
-  int nbOfCells=getNumberOfCells();
-  bool easyAssign=true;
-  const int *connI=_nodal_connec_index->getConstPointer();
-  const int *connIOther=otherOnSameCoordsThanThis._nodal_connec_index->getConstPointer();
+  std::size_t nbOfCells(getNumberOfCells());
+  bool easyAssign(true);
+  const int *connI(_nodal_connec_index->begin());
+  const int *connIOther=otherOnSameCoordsThanThis._nodal_connec_index->begin();
   for(const int *it=cellIdsBg;it!=cellIdsEnd && easyAssign;it++,connIOther++)
     {
-      if(*it>=0 && *it<nbOfCells)
+      if(*it>=0 && *it<(int)nbOfCells)
         {
           easyAssign=(connIOther[1]-connIOther[0])==(connI[*it+1]-connI[*it]);
         }
@@ -2034,7 +2042,7 @@ void MEDCouplingUMesh::setPartOfMySelfSlice(int start, int end, int step, const
       throw INTERP_KERNEL::Exception(oss.str());
     }
   int nbOfCellsToModify=DataArray::GetNumberOfItemGivenBESRelative(start,end,step,"MEDCouplingUMesh::setPartOfMySelfSlice : ");
-  if(nbOfCellsToModify!=otherOnSameCoordsThanThis.getNumberOfCells())
+  if(nbOfCellsToModify!=(int)otherOnSameCoordsThanThis.getNumberOfCells())
     {
       std::ostringstream oss; oss << "MEDCouplingUMesh::setPartOfMySelfSlice : cells ids length (" <<  nbOfCellsToModify << ") do not match the number of cells of other mesh (" << otherOnSameCoordsThanThis.getNumberOfCells() << ") !";
       throw INTERP_KERNEL::Exception(oss.str());
@@ -2277,7 +2285,7 @@ MEDCouplingUMesh *MEDCouplingUMesh::buildUnstructured() const
  * This method expects that \b this and \b otherDimM1OnSameCoords share the same coordinates array.
  * otherDimM1OnSameCoords->getMeshDimension() is expected to be equal to this->getMeshDimension()-1.
  * This method searches for nodes needed to be duplicated. These nodes are nodes fetched by \b otherDimM1OnSameCoords which are not part of the boundary of \b otherDimM1OnSameCoords.
- * If a node is in the boundary of \b this \b and in the boundary of \b otherDimM1OnSameCoords this node is considerd as needed to be duplicated.
+ * If a node is in the boundary of \b this \b and in the boundary of \b otherDimM1OnSameCoords this node is considered as needed to be duplicated.
  * When the set of node ids \b nodeIdsToDuplicate is computed, cell ids in \b this is searched so that their connectivity includes at least 1 node in \b nodeIdsToDuplicate.
  *
  * \param [in] otherDimM1OnSameCoords a mesh lying on the same coords than \b this and with a mesh dimension equal to those of \b this minus 1. WARNING this input
@@ -2363,7 +2371,7 @@ void MEDCouplingUMesh::findNodesToDuplicate(const MEDCouplingUMesh& otherDimM1On
   DAInt neighIInit00(tmp11);
   // Neighbor information of the mesh WITH the crack (some neighbors are removed):
   DataArrayInt *idsTmp=0;
-  bool b=m01->areCellsIncludedIn(&otherDimM1OnSameCoords,2,idsTmp);
+  m01->areCellsIncludedIn(&otherDimM1OnSameCoords,2,idsTmp);
   DAInt ids(idsTmp);
   // In the neighbor information remove the connection between high dimension cells and its low level constituents which are part
   // of the frontier given in parameter (i.e. the cells of low dimension from the group delimiting the crack):
@@ -2390,7 +2398,7 @@ void MEDCouplingUMesh::findNodesToDuplicate(const MEDCouplingUMesh& otherDimM1On
       // Connex zone without the crack (to compute the next seed really)
       int dnu;
       DAInt connexCheck = MEDCouplingUMesh::ComputeSpreadZoneGraduallyFromSeed(&seed, &seed+1, neighInit00,neighIInit00, -1, dnu);
-      int cnt = 0;
+      std::size_t cnt(0);
       for (int * ptr = connexCheck->getPointer(); cnt < connexCheck->getNumberOfTuples(); ptr++, cnt++)
         hitCells->setIJ(*ptr,0,1);
       // Connex zone WITH the crack (to identify cells lying on either part of the crack)
@@ -2422,7 +2430,7 @@ void MEDCouplingUMesh::findNodesToDuplicate(const MEDCouplingUMesh& otherDimM1On
  * This method operates a modification of the connectivity and coords in \b this.
  * Every time that a node id in [ \b nodeIdsToDuplicateBg, \b nodeIdsToDuplicateEnd ) will append in nodal connectivity of \b this 
  * its ids will be modified to id this->getNumberOfNodes()+std::distance(nodeIdsToDuplicateBg,std::find(nodeIdsToDuplicateBg,nodeIdsToDuplicateEnd,id)).
- * More explicitely the renumber array in nodes is not explicitely given in old2new to avoid to build a big array of renumbering whereas typically few node ids needs to be
+ * More explicitly the renumber array in nodes is not explicitly given in old2new to avoid to build a big array of renumbering whereas typically few node ids needs to be
  * renumbered. The node id nodeIdsToDuplicateBg[0] will have id this->getNumberOfNodes()+0, node id nodeIdsToDuplicateBg[1] will have id this->getNumberOfNodes()+1,
  * node id nodeIdsToDuplicateBg[2] will have id this->getNumberOfNodes()+2...
  * 
@@ -2563,7 +2571,7 @@ void MEDCouplingUMesh::shiftNodeNumbersInConn(int delta)
  * Coordinates are \b NOT considered here and will remain unchanged by this method. this->_coords can ever been null for the needs of this method.
  * Every time that a node id in [ \b nodeIdsToDuplicateBg, \b nodeIdsToDuplicateEnd ) will append in nodal connectivity of \b this 
  * its ids will be modified to id offset+std::distance(nodeIdsToDuplicateBg,std::find(nodeIdsToDuplicateBg,nodeIdsToDuplicateEnd,id)).
- * More explicitely the renumber array in nodes is not explicitely given in old2new to avoid to build a big array of renumbering whereas typically few node ids needs to be
+ * More explicitly the renumber array in nodes is not explicitly given in old2new to avoid to build a big array of renumbering whereas typically few node ids needs to be
  * renumbered. The node id nodeIdsToDuplicateBg[0] will have id offset+0, node id nodeIdsToDuplicateBg[1] will have id offset+1,
  * node id nodeIdsToDuplicateBg[2] will have id offset+2...
  * 
@@ -2775,11 +2783,10 @@ DataArrayInt *MEDCouplingUMesh::getCellsInBoundingBox(const INTERP_KERNEL::Direc
  *  \return INTERP_KERNEL::NormalizedCellType - enumeration item describing the cell type.
  *  \throw If \a cellId is invalid. Valid range is [0, \a this->getNumberOfCells() ).
  */
-INTERP_KERNEL::NormalizedCellType MEDCouplingUMesh::getTypeOfCell(int cellId) const
+INTERP_KERNEL::NormalizedCellType MEDCouplingUMesh::getTypeOfCell(std::size_t cellId) const
 {
-  const int *ptI=_nodal_connec_index->getConstPointer();
-  const int *pt=_nodal_connec->getConstPointer();
-  if(cellId>=0 && cellId<(int)_nodal_connec_index->getNbOfElems()-1)
+  const int *ptI(_nodal_connec_index->begin()),*pt(_nodal_connec->begin());
+  if(cellId<_nodal_connec_index->getNbOfElems()-1)
     return (INTERP_KERNEL::NormalizedCellType) pt[ptI[cellId]];
   else
     {
@@ -2821,13 +2828,11 @@ DataArrayInt *MEDCouplingUMesh::giveCellsWithType(INTERP_KERNEL::NormalizedCellT
 /*!
  * Returns nb of cells having the geometric type \a type. No throw if no cells in \a this has the geometric type \a type.
  */
-int MEDCouplingUMesh::getNumberOfCellsWithType(INTERP_KERNEL::NormalizedCellType type) const
+std::size_t MEDCouplingUMesh::getNumberOfCellsWithType(INTERP_KERNEL::NormalizedCellType type) const
 {
-  const int *ptI=_nodal_connec_index->getConstPointer();
-  const int *pt=_nodal_connec->getConstPointer();
-  int nbOfCells=getNumberOfCells();
-  int ret=0;
-  for(int i=0;i<nbOfCells;i++)
+  const int *ptI(_nodal_connec_index->begin()),*pt(_nodal_connec->begin());
+  std::size_t nbOfCells(getNumberOfCells()),ret(0);
+  for(std::size_t i=0;i<nbOfCells;i++)
     if((INTERP_KERNEL::NormalizedCellType) pt[ptI[i]]==type)
       ret++;
   return ret;
@@ -2842,10 +2847,9 @@ int MEDCouplingUMesh::getNumberOfCellsWithType(INTERP_KERNEL::NormalizedCellType
  *         cleared before the appending.
  *  \throw If \a cellId is invalid. Valid range is [0, \a this->getNumberOfCells() ).
  */
-void MEDCouplingUMesh::getNodeIdsOfCell(int cellId, std::vector<int>& conn) const
+void MEDCouplingUMesh::getNodeIdsOfCell(std::size_t cellId, std::vector<int>& conn) const
 {
-  const int *ptI=_nodal_connec_index->getConstPointer();
-  const int *pt=_nodal_connec->getConstPointer();
+  const int *ptI(_nodal_connec_index->begin()),*pt(_nodal_connec->begin());
   for(const int *w=pt+ptI[cellId]+1;w!=pt+ptI[cellId+1];w++)
     if(*w>=0)
       conn.push_back(*w);
@@ -2938,7 +2942,7 @@ std::string MEDCouplingUMesh::reprConnectivityOfThis() const
 }
 
 /*!
- * This method builds a newly allocated instance (with the same name than \a this) that the caller has the responsability to deal with.
+ * This method builds a newly allocated instance (with the same name than \a this) that the caller has the responsibility to deal with.
  * This method returns an instance with all arrays allocated (connectivity, connectivity index, coordinates)
  * but with length of these arrays set to 0. It allows to define an "empty" mesh (with nor cells nor nodes but compliant with
  * some algos).
@@ -3000,7 +3004,7 @@ int MEDCouplingUMesh::getNumberOfNodesInCell(int cellId) const
 
 /*!
  * Returns types of cells of the specified part of \a this mesh.
- * This method avoids computing sub-mesh explicitely to get its types.
+ * This method avoids computing sub-mesh explicitly to get its types.
  *  \param [in] begin - an array of cell ids of interest.
  *  \param [in] end - the end of \a begin, i.e. a pointer to its (last+1)-th element.
  *  \return std::set<INTERP_KERNEL::NormalizedCellType> - a set of enumeration items
@@ -3044,14 +3048,14 @@ void MEDCouplingUMesh::setConnectivity(DataArrayInt *conn, DataArrayInt *connInd
  * Copy constructor. If 'deepCopy' is false \a this is a shallow copy of other.
  * If 'deeCpy' is true all arrays (coordinates and connectivities) are deeply copied.
  */
-MEDCouplingUMesh::MEDCouplingUMesh(const MEDCouplingUMesh& other, bool deepCopy):MEDCouplingPointSet(other,deepCopy),_mesh_dim(other._mesh_dim),
+MEDCouplingUMesh::MEDCouplingUMesh(const MEDCouplingUMesh& other, bool deepCpy):MEDCouplingPointSet(other,deepCpy),_mesh_dim(other._mesh_dim),
     _nodal_connec(0),_nodal_connec_index(0),
     _types(other._types)
 {
   if(other._nodal_connec)
-    _nodal_connec=other._nodal_connec->performCopyOrIncrRef(deepCopy);
+    _nodal_connec=other._nodal_connec->performCopyOrIncrRef(deepCpy);
   if(other._nodal_connec_index)
-    _nodal_connec_index=other._nodal_connec_index->performCopyOrIncrRef(deepCopy);
+    _nodal_connec_index=other._nodal_connec_index->performCopyOrIncrRef(deepCpy);
 }
 
 MEDCouplingUMesh::~MEDCouplingUMesh()
@@ -3843,7 +3847,7 @@ DataArrayInt *MEDCouplingUMesh::getCellIdsCrossingPlane(const double *origin, co
  * If not an exception will thrown. If this is an empty mesh with no cell an exception will be thrown too.
  * No consideration of coordinate is done by this method.
  * 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 MEDCoupling::MEDCouplingUMesh::mergeNodes could be usefull.
+ * If not false is returned. In case that false is returned a call to MEDCoupling::MEDCouplingUMesh::mergeNodes could be useful.
  */
 bool MEDCouplingUMesh::isContiguous1D() const
 {
@@ -3882,7 +3886,7 @@ void MEDCouplingUMesh::project1D(const double *pt, const double *v, double eps,
   MCAuto<MEDCouplingFieldDouble> f=buildDirectionVectorField();
   const double *fPtr=f->getArray()->getConstPointer();
   double tmp[3];
-  for(int i=0;i<getNumberOfCells();i++)
+  for(std::size_t i=0;i<getNumberOfCells();i++)
     {
       const double *tmp1=fPtr+3*i;
       tmp[0]=tmp1[1]*v[2]-tmp1[2]*v[1];
@@ -3972,7 +3976,7 @@ DataArrayDouble *MEDCouplingUMesh::distanceToPoints(const DataArrayDouble *pts,
     throw INTERP_KERNEL::Exception("MEDCouplingUMesh::distanceToPoints works only for spaceDim=meshDim+1 !");
   if(meshDim!=2 && meshDim!=1)
     throw INTERP_KERNEL::Exception("MEDCouplingUMesh::distanceToPoints : only mesh dimension 2 and 1 are implemented !");
-  if(pts->getNumberOfComponents()!=spaceDim)
+  if((int)pts->getNumberOfComponents()!=spaceDim)
     {
       std::ostringstream oss; oss << "MEDCouplingUMesh::distanceToPoints : input pts DataArrayDouble has " << pts->getNumberOfComponents() << " components whereas it should be equal to " << spaceDim << " (mesh spaceDimension) !";
       throw INTERP_KERNEL::Exception(oss.str());
@@ -5421,6 +5425,8 @@ DataArrayDouble *MEDCouplingUMesh::getBoundingBoxForBBTreeFast() const
 DataArrayDouble *MEDCouplingUMesh::getBoundingBoxForBBTree2DQuadratic(double arcDetEps) const
 {
   checkFullyDefined();
+  INTERP_KERNEL::QuadraticPlanarArcDetectionPrecision arcPrec(arcDetEps);
+
   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!");
@@ -5432,7 +5438,6 @@ DataArrayDouble *MEDCouplingUMesh::getBoundingBoxForBBTree2DQuadratic(double arc
     {
       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++)
@@ -5469,6 +5474,7 @@ DataArrayDouble *MEDCouplingUMesh::getBoundingBoxForBBTree1DQuadratic(double arc
   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!");
+  INTERP_KERNEL::QuadraticPlanarArcDetectionPrecision arcPrec(arcDetEps);
   MCAuto<DataArrayDouble> ret(DataArrayDouble::New()); ret->alloc(nbOfCells,2*spaceDim);
   double *bbox(ret->getPointer());
   const double *coords(_coords->begin());
@@ -5477,7 +5483,6 @@ DataArrayDouble *MEDCouplingUMesh::getBoundingBoxForBBTree1DQuadratic(double arc
     {
       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++)
@@ -5528,7 +5533,7 @@ namespace MEDCouplingImpl
  * \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]==-1 because it has no sense here. 
- * This parameter is kept only for compatibility with other methode listed above.
+ * This parameter is kept only for compatibility with other method listed above.
  */
 std::vector<int> MEDCouplingUMesh::getDistributionOfTypes() const
 {
@@ -5573,7 +5578,7 @@ std::vector<int> MEDCouplingUMesh::getDistributionOfTypes() const
  * 
  * If all geometric types in \a code are exactly those in \a this null pointer is returned.
  * 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.
+ * and a DataArrayInt instance is returned that the user has the responsibility to deallocate.
  */
 DataArrayInt *MEDCouplingUMesh::checkTypeConsistencyAndContig(const std::vector<int>& code, const std::vector<const DataArrayInt *>& idsPerType) const
 {
@@ -5884,7 +5889,7 @@ bool MEDCouplingUMesh::checkConsecutiveCellTypesAndOrder(const INTERP_KERNEL::No
 /*!
  * This method returns 2 newly allocated DataArrayInt instances. The first is an array of size 'this->getNumberOfCells()' with one component,
  * 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'.
+ * number of tuples than input type array and with one component. This 2nd output array gives type by type the number of occurrence of type in 'this'.
  */
 DataArrayInt *MEDCouplingUMesh::getLevArrPerCellTypes(const INTERP_KERNEL::NormalizedCellType *orderBg, const INTERP_KERNEL::NormalizedCellType *orderEnd, DataArrayInt *&nbPerType) const
 {
@@ -6125,7 +6130,7 @@ void MEDCouplingUMesh::convertNodalConnectivityToDynamicGeoTypeMesh(DataArrayInt
 /*!
  * 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).
- * This method is particulary useful for MED file interaction. It allows to aggregate several meshes and keeping the type sorting
+ * This method is particularly useful for MED file interaction. It allows to aggregate several meshes and keeping the type sorting
  * and the track of the permutation by chunk of same geotype cells to retrieve it. The traditional formats old2new and new2old
  * are not used here to avoid the build of big permutation array.
  *
@@ -6479,8 +6484,15 @@ DataArrayDouble *MEDCouplingUMesh::computePlaneEquationOf3DFaces() const
   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)
+      if(nodalI[1]-nodalI[0]>=4)
         {
+          double aa[3]={coor[nodal[nodalI[0]+1+1]*3+0]-coor[nodal[nodalI[0]+1+0]*3+0],
+                        coor[nodal[nodalI[0]+1+1]*3+1]-coor[nodal[nodalI[0]+1+0]*3+1],
+                        coor[nodal[nodalI[0]+1+1]*3+2]-coor[nodal[nodalI[0]+1+0]*3+2]}
+          ,bb[3]={coor[nodal[nodalI[0]+1+2]*3+0]-coor[nodal[nodalI[0]+1+0]*3+0],
+                        coor[nodal[nodalI[0]+1+2]*3+1]-coor[nodal[nodalI[0]+1+0]*3+1],
+                        coor[nodal[nodalI[0]+1+2]*3+2]-coor[nodal[nodalI[0]+1+0]*3+2]};
+          double cc[3]={aa[1]*bb[2]-aa[2]*bb[1],aa[2]*bb[0]-aa[0]*bb[2],aa[0]*bb[1]-aa[1]*bb[0]};
           for(int j=0;j<3;j++)
             {
               int nodeId(nodal[nodalI[0]+1+j]);
@@ -6492,14 +6504,34 @@ DataArrayDouble *MEDCouplingUMesh::computePlaneEquationOf3DFaces() const
                   throw INTERP_KERNEL::Exception(oss.str());
                 }
             }
+          if(sqrt(cc[0]*cc[0]+cc[1]*cc[1]+cc[2]*cc[2])>1e-7)
+            {
+              INTERP_KERNEL::inverseMatrix(matrix,4,matrix2);
+              retPtr[0]=matrix2[3]; retPtr[1]=matrix2[7]; retPtr[2]=matrix2[11]; retPtr[3]=matrix2[15];
+            }
+          else
+            {
+              if(nodalI[1]-nodalI[0]==4)
+                {
+                  std::ostringstream oss; oss << "MEDCouplingUMesh::computePlaneEquationOf3DFaces : cell" << i << " : Presence of The 3 colinear points !";
+                  throw INTERP_KERNEL::Exception(oss.str());
+                }
+              //
+              double dd[3]={0.,0.,0.};
+              for(int offset=nodalI[0]+1;offset<nodalI[1];offset++)
+                std::transform(coor+3*nodal[offset],coor+3*(nodal[offset]+1),dd,dd,std::plus<double>());
+              int nbOfNodesInCell(nodalI[1]-nodalI[0]-1);
+              std::transform(dd,dd+3,dd,std::bind2nd(std::multiplies<double>(),1./(double)nbOfNodesInCell));
+              std::copy(dd,dd+3,matrix+4*2);
+              INTERP_KERNEL::inverseMatrix(matrix,4,matrix2);
+              retPtr[0]=matrix2[3]; retPtr[1]=matrix2[7]; retPtr[2]=matrix2[11]; retPtr[3]=matrix2[15];
+            }
         }
       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());
         }
-      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();
 }
@@ -6761,7 +6793,7 @@ MEDCouplingUMesh *MEDCouplingUMesh::FuseUMeshesOnSameCoords(const std::vector<co
  * Makes all given meshes share the nodal connectivity array. The common connectivity
  * array is created by concatenating the connectivity arrays of all given meshes. All
  * the given meshes must be of the same space dimension but dimension of cells **can
- * differ**. This method is particulary useful in MEDLoader context to build a \ref
+ * differ**. This method is particularly useful in MEDLoader context to build a \ref
  * MEDCoupling::MEDFileUMesh "MEDFileUMesh" instance that expects that underlying
  * MEDCouplingUMesh'es of different dimension share the same nodal connectivity array.
  *  \param [in,out] meshes - a vector of meshes to update.
@@ -6815,7 +6847,7 @@ void MEDCouplingUMesh::PutUMeshesOnSameAggregatedCoords(const std::vector<MEDCou
 /*!
  * Merges nodes coincident with a given precision within all given meshes that share
  * the nodal connectivity array. The given meshes **can be of different** mesh
- * dimension. This method is particulary useful in MEDLoader context to build a \ref
+ * dimension. This method is particularly useful in MEDLoader context to build a \ref
  * MEDCoupling::MEDFileUMesh "MEDFileUMesh" instance that expects that underlying
  * MEDCouplingUMesh'es of different dimension share the same nodal connectivity array. 
  *  \param [in,out] meshes - a vector of meshes to update.
@@ -6996,48 +7028,62 @@ 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)
+void MEDCouplingUMesh::SimplifyPolyhedronCell(double eps, const DataArrayDouble *coords, int index, DataArrayInt *res, MEDCouplingUMesh *faces,
+                                              DataArrayInt *E_Fi, DataArrayInt *E_F, DataArrayInt *F_Ei, DataArrayInt *F_E)
 {
-  int nbFaces=std::count(begin+1,end,-1)+1;
+  int nbFaces = E_Fi->getIJ(index + 1, 0) - E_Fi->getIJ(index, 0);
   MCAuto<DataArrayDouble> v=DataArrayDouble::New(); v->alloc(nbFaces,3);
   double *vPtr=v->getPointer();
-  MCAuto<DataArrayDouble> p=DataArrayDouble::New(); p->alloc(nbFaces,1);
+  MCAuto<DataArrayDouble> p=DataArrayDouble::New(); p->alloc(nbFaces,2);
   double *pPtr=p->getPointer();
-  const int *stFaceConn=begin+1;
+  int *e_fi = E_Fi->getPointer(), *e_f = E_F->getPointer(), *f_ei = F_Ei->getPointer(), *f_e = F_E->getPointer();
+  const int *f_idx = faces->getNodalConnectivityIndex()->getPointer(), *f_cnn = faces->getNodalConnectivity()->getPointer();
   for(int i=0;i<nbFaces;i++,vPtr+=3,pPtr++)
     {
-      const int *endFaceConn=std::find(stFaceConn,end,-1);
-      ComputeVecAndPtOfFace(eps,coords->begin(),stFaceConn,endFaceConn,vPtr,pPtr);
-      stFaceConn=endFaceConn+1;
+      int face = e_f[e_fi[index] + i];
+      ComputeVecAndPtOfFace(eps, coords->begin(), f_cnn + f_idx[face] + 1, f_cnn + f_idx[face + 1], vPtr, pPtr);
+      // to differentiate faces going to different cells:
+      pPtr++, *pPtr = 0;
+      for (int j = f_ei[face]; j < f_ei[face + 1]; j++)
+        *pPtr += f_e[j];
     }
   pPtr=p->getPointer(); vPtr=v->getPointer();
   DataArrayInt *comm1=0,*commI1=0;
   v->findCommonTuples(eps,-1,comm1,commI1);
+  for (int i = 0; i < nbFaces; i++)
+    if (comm1->findIdFirstEqual(i) < 0)
+      {
+        comm1->pushBackSilent(i);
+        commI1->pushBackSilent(comm1->getNumberOfTuples());
+      }
   MCAuto<DataArrayInt> comm1Auto(comm1),commI1Auto(commI1);
   const int *comm1Ptr=comm1->begin();
   const int *commI1Ptr=commI1->begin();
   int nbOfGrps1=commI1Auto->getNumberOfTuples()-1;
   res->pushBackSilent((int)INTERP_KERNEL::NORM_POLYHED);
   //
-  MCAuto<MEDCouplingUMesh> mm=MEDCouplingUMesh::New("",3);
-  mm->setCoords(const_cast<DataArrayDouble *>(coords)); mm->allocateCells(1); mm->insertNextCell(INTERP_KERNEL::NORM_POLYHED,(int)std::distance(begin+1,end),begin+1);
-  mm->finishInsertingCells();
-  //
   for(int i=0;i<nbOfGrps1;i++)
     {
       int vecId=comm1Ptr[commI1Ptr[i]];
       MCAuto<DataArrayDouble> tmpgrp2=p->selectByTupleId(comm1Ptr+commI1Ptr[i],comm1Ptr+commI1Ptr[i+1]);
       DataArrayInt *comm2=0,*commI2=0;
       tmpgrp2->findCommonTuples(eps,-1,comm2,commI2);
+      for (int j = 0; j < commI1Ptr[i+1] - commI1Ptr[i]; j++)
+        if (comm2->findIdFirstEqual(j) < 0)
+          {
+            comm2->pushBackSilent(j);
+            commI2->pushBackSilent(comm2->getNumberOfTuples());
+          }
       MCAuto<DataArrayInt> comm2Auto(comm2),commI2Auto(commI2);
       const int *comm2Ptr=comm2->begin();
       const int *commI2Ptr=commI2->begin();
       int nbOfGrps2=commI2Auto->getNumberOfTuples()-1;
       for(int j=0;j<nbOfGrps2;j++)
         {
-          if(commI2Ptr[j+1]-commI2Ptr[j]<=1)
+          if(commI2Ptr[j+1] == commI2Ptr[j] + 1)
             {
-              res->insertAtTheEnd(begin,end);
+              int face = e_f[e_fi[index] + comm1Ptr[commI1Ptr[i] + comm2Ptr[commI2Ptr[j]]]]; //hmmm
+              res->insertAtTheEnd(f_cnn + f_idx[face] + 1, f_cnn + f_idx[face + 1]);
               res->pushBackSilent(-1);
             }
           else
@@ -7045,13 +7091,12 @@ void MEDCouplingUMesh::SimplifyPolyhedronCell(double eps, const DataArrayDouble
               int pointId=comm1Ptr[commI1Ptr[i]+comm2Ptr[commI2Ptr[j]]];
               MCAuto<DataArrayInt> ids2=comm2->selectByTupleIdSafeSlice(commI2Ptr[j],commI2Ptr[j+1],1);
               ids2->transformWithIndArr(comm1Ptr+commI1Ptr[i],comm1Ptr+commI1Ptr[i+1]);
-              DataArrayInt *tmp0=DataArrayInt::New(),*tmp1=DataArrayInt::New(),*tmp2=DataArrayInt::New(),*tmp3=DataArrayInt::New();
-              MCAuto<MEDCouplingUMesh> mm2=mm->buildDescendingConnectivity(tmp0,tmp1,tmp2,tmp3); tmp0->decrRef(); tmp1->decrRef(); tmp2->decrRef(); tmp3->decrRef();
-              MCAuto<MEDCouplingUMesh> mm3=static_cast<MEDCouplingUMesh *>(mm2->buildPartOfMySelf(ids2->begin(),ids2->end(),true));
+              ids2->transformWithIndArr(e_f + e_fi[index], e_f + e_fi[index + 1]);
+              MCAuto<MEDCouplingUMesh> mm3=static_cast<MEDCouplingUMesh *>(faces->buildPartOfMySelf(ids2->begin(),ids2->end(),true));
               MCAuto<DataArrayInt> idsNodeTmp=mm3->zipCoordsTraducer();
               MCAuto<DataArrayInt> idsNode=idsNodeTmp->invertArrayO2N2N2O(mm3->getNumberOfNodes());
               const int *idsNodePtr=idsNode->begin();
-              double center[3]; center[0]=pPtr[pointId]*vPtr[3*vecId]; center[1]=pPtr[pointId]*vPtr[3*vecId+1]; center[2]=pPtr[pointId]*vPtr[3*vecId+2];
+              double center[3]; center[0]=pPtr[2*pointId]*vPtr[3*vecId]; center[1]=pPtr[2*pointId]*vPtr[3*vecId+1]; center[2]=pPtr[2*pointId]*vPtr[3*vecId+2];
               double vec[3]; vec[0]=vPtr[3*vecId+1]; vec[1]=-vPtr[3*vecId]; vec[2]=0.;
               double norm=vec[0]*vec[0]+vec[1]*vec[1]+vec[2]*vec[2];
               if(std::abs(norm)>eps)
@@ -7703,7 +7748,7 @@ bool MEDCouplingUMesh::RemoveIdsFromIndexedArrays(const int *idsToRemoveBg, cons
       previousArrI=*arrIPtr;
       *arrIPtr=(int)arrOut.size();
     }
-  if(arr->getNumberOfTuples()==(int)arrOut.size())
+  if(arr->getNumberOfTuples()==arrOut.size())
     return false;
   arr->alloc((int)arrOut.size(),1);
   std::copy(arrOut.begin(),arrOut.end(),arr->getPointer());
@@ -7861,7 +7906,7 @@ void MEDCouplingUMesh::ExtractFromIndexedArraysSlice(int idsOfSelectStart, int i
  * This method works on an input pair (\b arrIn, \b arrIndxIn) where \b arrIn indexes is in \b arrIndxIn.
  * This method builds an output pair (\b arrOut,\b arrIndexOut) that is a copy from \b arrIn for all cell ids \b not \b in [ \b idsOfSelectBg , \b idsOfSelectEnd ) and for
  * cellIds \b in [ \b idsOfSelectBg , \b idsOfSelectEnd ) a copy coming from the corresponding values in input pair (\b srcArr, \b srcArrIndex).
- * This method is an generalization of MEDCouplingUMesh::SetPartOfIndexedArraysSameIdx that performs the same thing but by without building explicitely a result output arrays.
+ * This method is an generalization of MEDCouplingUMesh::SetPartOfIndexedArraysSameIdx that performs the same thing but by without building explicitly a result output arrays.
  *
  * \param [in] idsOfSelectBg begin of set of ids of the input extraction (included)
  * \param [in] idsOfSelectEnd end of set of ids of the input extraction (excluded)
@@ -7927,7 +7972,7 @@ void MEDCouplingUMesh::SetPartOfIndexedArrays(const int *idsOfSelectBg, const in
 
 /*!
  * This method works on an input pair (\b arrIn, \b arrIndxIn) where \b arrIn indexes is in \b arrIndxIn.
- * This method is an specialization of MEDCouplingUMesh::SetPartOfIndexedArrays in the case of assignement do not modify the index in \b arrIndxIn.
+ * This method is an specialization of MEDCouplingUMesh::SetPartOfIndexedArrays in the case of assignment do not modify the index in \b arrIndxIn.
  *
  * \param [in] idsOfSelectBg begin of set of ids of the input extraction (included)
  * \param [in] idsOfSelectEnd end of set of ids of the input extraction (excluded)
@@ -8024,7 +8069,7 @@ DataArrayInt *MEDCouplingUMesh::ComputeSpreadZoneGraduallyFromSeed(const int *se
  * This method works on an input pair (\b arrIn, \b arrIndxIn) where \b arrIn indexes is in \b arrIndxIn.
  * This method builds an output pair (\b arrOut,\b arrIndexOut) that is a copy from \b arrIn for all cell ids \b not \b in [ \b idsOfSelectBg , \b idsOfSelectEnd ) and for
  * cellIds \b in [\b idsOfSelectBg, \b idsOfSelectEnd) a copy coming from the corresponding values in input pair (\b srcArr, \b srcArrIndex).
- * This method is an generalization of MEDCouplingUMesh::SetPartOfIndexedArraysSameIdx that performs the same thing but by without building explicitely a result output arrays.
+ * This method is an generalization of MEDCouplingUMesh::SetPartOfIndexedArraysSameIdx that performs the same thing but by without building explicitly a result output arrays.
  *
  * \param [in] start begin of set of ids of the input extraction (included)
  * \param [in] end end of set of ids of the input extraction (excluded)
@@ -8089,7 +8134,7 @@ void MEDCouplingUMesh::SetPartOfIndexedArraysSlice(int start, int end, int step,
 
 /*!
  * This method works on an input pair (\b arrIn, \b arrIndxIn) where \b arrIn indexes is in \b arrIndxIn.
- * This method is an specialization of MEDCouplingUMesh::SetPartOfIndexedArrays in the case of assignement do not modify the index in \b arrIndxIn.
+ * This method is an specialization of MEDCouplingUMesh::SetPartOfIndexedArrays in the case of assignment do not modify the index in \b arrIndxIn.
  *
  * \param [in] start begin of set of ids of the input extraction (included)
  * \param [in] end end of set of ids of the input extraction (excluded)
@@ -8259,7 +8304,7 @@ DataArrayInt *MEDCouplingUMesh::ComputeRangesFromTypeDistribution(const std::vec
  *         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
+ * \warning This method operates on each cells in this independently ! 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).