Salome HOME
Copyright update 2022
[tools/medcoupling.git] / src / INTERP_KERNEL / Geometric2D / InterpKernelGeo2DQuadraticPolygon.cxx
index 9ada88d0406737c1166a4c5aeeced6dbb4098525..f6c1ee3164a5cbcf525d754ca07fad61c474a9b4 100644 (file)
@@ -1,4 +1,4 @@
-// Copyright (C) 2007-2016  CEA/DEN, EDF R&D
+// Copyright (C) 2007-2022  CEA/DEN, EDF R&D
 //
 // This library is free software; you can redistribute it and/or
 // modify it under the terms of the GNU Lesser General Public
@@ -284,12 +284,12 @@ double QuadraticPolygon::intersectWithAbs(QuadraticPolygon& other)
  * the cell id in global other mesh.
  */
 void QuadraticPolygon::splitAbs(QuadraticPolygon& other,
-                                const std::map<INTERP_KERNEL::Node *,int>& mapThis, const std::map<INTERP_KERNEL::Node *,int>& mapOther,
-                                int offset1, int offset2 ,
-                                const std::vector<int>& otherEdgeIds,
-                                std::vector<int>& edgesThis, int cellIdThis,
-                                std::vector< std::vector<int> >& edgesInOtherColinearWithThis, std::vector< std::vector<int> >& subDivOther,
-                                std::vector<double>& addCoo, std::map<int,int>& mergedNodes)
+                                const std::map<INTERP_KERNEL::Node *,mcIdType>& mapThis, const std::map<INTERP_KERNEL::Node *,mcIdType>& mapOther,
+                                mcIdType offset1, mcIdType offset2 ,
+                                const std::vector<mcIdType>& otherEdgeIds,
+                                std::vector<mcIdType>& edgesThis, mcIdType cellIdThis,
+                                std::vector< std::vector<mcIdType> >& edgesInOtherColinearWithThis, std::vector< std::vector<mcIdType> >& subDivOther,
+                                std::vector<double>& addCoo, std::map<mcIdType,mcIdType>& mergedNodes)
 {
   double xBaryBB, yBaryBB;
   double fact=normalizeExt(&other, xBaryBB, yBaryBB);
@@ -300,7 +300,7 @@ void QuadraticPolygon::splitAbs(QuadraticPolygon& other,
   ComposedEdge *cThis=new ComposedEdge;
   ComposedEdge *cOther=new ComposedEdge;
   int i=0;
-  std::map<INTERP_KERNEL::Node *,int> mapAddCoo;
+  std::map<INTERP_KERNEL::Node *,mcIdType> mapAddCoo;
   for(itOther.first();!itOther.finished();itOther.next(),i++)
     {
       // For each edge of 'other', proceed with intersections: the edge might split into sub-edges, 'otherTmp' will hold the final split result.
@@ -322,9 +322,9 @@ void QuadraticPolygon::splitAbs(QuadraticPolygon& other,
               ElementaryEdge* curThis=itThis.current();
               merge.clear();
               //
-              std::map<INTERP_KERNEL::Node *,int>::const_iterator thisStart(mapThis.find(curThis->getStartNode())),thisEnd(mapThis.find(curThis->getEndNode())),
+              std::map<INTERP_KERNEL::Node *,mcIdType>::const_iterator thisStart(mapThis.find(curThis->getStartNode())),thisEnd(mapThis.find(curThis->getEndNode())),
                                                                   otherStart(mapOther.find(curOtherTmp->getStartNode())),otherEnd(mapOther.find(curOtherTmp->getEndNode()));
-              int thisStart2(thisStart==mapThis.end()?-1:(*thisStart).second), thisEnd2(thisEnd==mapThis.end()?-1:(*thisEnd).second),
+              mcIdType thisStart2(thisStart==mapThis.end()?-1:(*thisStart).second), thisEnd2(thisEnd==mapThis.end()?-1:(*thisEnd).second),
                   otherStart2(otherStart==mapOther.end()?-1:(*otherStart).second+offset1),otherEnd2(otherEnd==mapOther.end()?-1:(*otherEnd).second+offset1);
               //
               if(curThis->getPtr()->intersectWith(curOtherTmp->getPtr(),merge,*cThis,*cOther))
@@ -362,10 +362,10 @@ void QuadraticPolygon::splitAbs(QuadraticPolygon& other,
       // Converting back to integer connectivity:
       if(otherTmp._sub_edges.size()>1)   // only if a new point has been added (i.e. an actual intersection was done)
         {
-          int jj = 0;
+          std::size_t jj = 0, sz(otherTmp._sub_edges.size());
           for(std::list<ElementaryEdge *>::const_iterator it=otherTmp._sub_edges.begin();it!=otherTmp._sub_edges.end();it++, jj++)
             {
-              unsigned skipStartOrEnd = jj == 0 ? 1 : (jj == _sub_edges.size()-1 ? 2 : -1);  // 1 means START, 2 means END, -1 other
+              short skipStartOrEnd = jj == 0 ? -1 : (jj == sz-1 ? 1 : 0);  // -1 means START, 1 means END, 0 other
               (*it)->fillGlobalInfoAbs2(mapThis,mapOther,offset1,offset2,
                                       fact,xBaryBB,yBaryBB, skipStartOrEnd,
                                       /*out*/ subDivOther[otherEdgeIds[i]],addCoo,mapAddCoo);
@@ -386,8 +386,8 @@ void QuadraticPolygon::splitAbs(QuadraticPolygon& other,
  * orientation of edge (see buildDescendingConnectivity2() method).
  * See appendEdgeFromCrudeDataArray() for params description.
  */
-void QuadraticPolygon::buildFromCrudeDataArray(const std::map<int,INTERP_KERNEL::Node *>& mapp, bool isQuad, const int *nodalBg, const double *coords,
-                                               const int *descBg, const int *descEnd, const std::vector<std::vector<int> >& intersectEdges)
+void QuadraticPolygon::buildFromCrudeDataArray(const std::map<mcIdType,INTERP_KERNEL::Node *>& mapp, bool isQuad, const mcIdType *nodalBg, const double *coords,
+                                               const mcIdType *descBg, const mcIdType *descEnd, const std::vector<std::vector<mcIdType> >& intersectEdges)
 {
   std::size_t nbOfSeg=std::distance(descBg,descEnd);
   for(std::size_t i=0;i<nbOfSeg;i++)
@@ -396,15 +396,15 @@ void QuadraticPolygon::buildFromCrudeDataArray(const std::map<int,INTERP_KERNEL:
     }
 }
 
-void QuadraticPolygon::appendEdgeFromCrudeDataArray(std::size_t edgePos, const std::map<int,INTERP_KERNEL::Node *>& mapp, bool isQuad,
-                                                    const int *nodalBg, const double *coords,
-                                                    const int *descBg, const int *descEnd, const std::vector<std::vector<int> >& intersectEdges)
+void QuadraticPolygon::appendEdgeFromCrudeDataArray(std::size_t edgePos, const std::map<mcIdType,INTERP_KERNEL::Node *>& mapp, bool isQuad,
+                                                    const mcIdType *nodalBg, const double *coords,
+                                                    const mcIdType *descBg, const mcIdType *descEnd, const std::vector<std::vector<mcIdType> >& intersectEdges)
 {
   if(!isQuad)
     {
       bool direct=descBg[edgePos]>0;
-      int edgeId=abs(descBg[edgePos])-1; // back to C indexing mode
-      const std::vector<int>& subEdge=intersectEdges[edgeId];
+      mcIdType edgeId=std::abs(descBg[edgePos])-1; // back to C indexing mode
+      const std::vector<mcIdType>& subEdge=intersectEdges[edgeId];
       std::size_t nbOfSubEdges=subEdge.size()/2;
       for(std::size_t j=0;j<nbOfSubEdges;j++)
         appendSubEdgeFromCrudeDataArray(0,j,direct,edgeId,subEdge,mapp);
@@ -426,8 +426,8 @@ void QuadraticPolygon::appendEdgeFromCrudeDataArray(std::size_t edgePos, const s
       delete e1; delete e2;
       //
       bool direct=descBg[edgePos]>0;
-      int edgeId=abs(descBg[edgePos])-1;
-      const std::vector<int>& subEdge=intersectEdges[edgeId];
+      mcIdType edgeId=std::abs(descBg[edgePos])-1;
+      const std::vector<mcIdType>& subEdge=intersectEdges[edgeId];
       std::size_t nbOfSubEdges=subEdge.size()/2;
       if(colinearity)
         {   
@@ -445,7 +445,7 @@ void QuadraticPolygon::appendEdgeFromCrudeDataArray(std::size_t edgePos, const s
     }
 }
 
-void QuadraticPolygon::appendSubEdgeFromCrudeDataArray(Edge *baseEdge, std::size_t j, bool direct, int edgeId, const std::vector<int>& subEdge, const std::map<int,INTERP_KERNEL::Node *>& mapp)
+void QuadraticPolygon::appendSubEdgeFromCrudeDataArray(Edge *baseEdge, std::size_t j, bool direct, mcIdType edgeId, const std::vector<mcIdType>& subEdge, const std::map<mcIdType,INTERP_KERNEL::Node *>& mapp)
 {
   std::size_t nbOfSubEdges=subEdge.size()/2;
   if(!baseEdge)
@@ -469,17 +469,17 @@ void QuadraticPolygon::appendSubEdgeFromCrudeDataArray(Edge *baseEdge, std::size
  * This method builds from descending conn of a quadratic polygon stored in crude mode (MEDCoupling). Descending conn is in FORTRAN relative mode in order to give the
  * orientation of edge.
  */
-void QuadraticPolygon::buildFromCrudeDataArray2(const std::map<int,INTERP_KERNEL::Node *>& mapp, bool isQuad, const int *nodalBg, const double *coords, const int *descBg, const int *descEnd, const std::vector<std::vector<int> >& intersectEdges2,
-                                                const INTERP_KERNEL::QuadraticPolygon& pol1, const int *descBg1, const int *descEnd1, const std::vector<std::vector<int> >& intersectEdges1,
-                                                const std::vector< std::vector<int> >& colinear1,
-                                                std::map<int,std::vector<INTERP_KERNEL::ElementaryEdge *> >& alreadyExistingIn2)
+void QuadraticPolygon::buildFromCrudeDataArray2(const std::map<mcIdType,INTERP_KERNEL::Node *>& mapp, bool isQuad, const mcIdType *nodalBg, const double *coords, const mcIdType *descBg, const mcIdType *descEnd, const std::vector<std::vector<mcIdType> >& intersectEdges2,
+                                                const INTERP_KERNEL::QuadraticPolygon& pol1, const mcIdType *descBg1, const mcIdType *descEnd1, const std::vector<std::vector<mcIdType> >& intersectEdges1,
+                                                const std::vector< std::vector<mcIdType> >& colinear1,
+                                                std::map<mcIdType,std::vector<INTERP_KERNEL::ElementaryEdge *> >& alreadyExistingIn2)
 {
   std::size_t nbOfSeg=std::distance(descBg,descEnd);
   for(std::size_t i=0;i<nbOfSeg;i++)//loop over all edges of pol2
     {
       bool direct=descBg[i]>0;
-      int edgeId=abs(descBg[i])-1;//current edge id of pol2
-      std::map<int,std::vector<INTERP_KERNEL::ElementaryEdge *> >::const_iterator it1=alreadyExistingIn2.find(descBg[i]),it2=alreadyExistingIn2.find(-descBg[i]);
+      mcIdType edgeId=std::abs(descBg[i])-1;//current edge id of pol2
+      std::map<mcIdType,std::vector<INTERP_KERNEL::ElementaryEdge *> >::const_iterator it1=alreadyExistingIn2.find(descBg[i]),it2=alreadyExistingIn2.find(-descBg[i]);
       if(it1!=alreadyExistingIn2.end() || it2!=alreadyExistingIn2.end())
         {
           bool sameDir=(it1!=alreadyExistingIn2.end());
@@ -503,21 +503,21 @@ void QuadraticPolygon::buildFromCrudeDataArray2(const std::map<int,INTERP_KERNEL
           continue;
         }
       bool directos=colinear1[edgeId].empty();
-      std::vector<std::pair<int,std::pair<bool,int> > > idIns1;
-      int offset1=0;
+      std::vector<std::pair<mcIdType,std::pair<bool,mcIdType> > > idIns1;
+      mcIdType offset1=0;
       if(!directos)
         {// if the current edge of pol2 has one or more colinear edges part into pol1
-          const std::vector<int>& c=colinear1[edgeId];
+          const std::vector<mcIdType>& c=colinear1[edgeId];
           std::size_t nbOfEdgesIn1=std::distance(descBg1,descEnd1);
           for(std::size_t j=0;j<nbOfEdgesIn1;j++)
             {
-              int edgeId1=abs(descBg1[j])-1;
+              mcIdType edgeId1=std::abs(descBg1[j])-1;
               if(std::find(c.begin(),c.end(),edgeId1)!=c.end())
                 {
-                  idIns1.push_back(std::pair<int,std::pair<bool,int> >(edgeId1,std::pair<bool,int>(descBg1[j]>0,offset1)));// it exists an edge into pol1 given by tuple (idIn1,direct1) that is colinear at edge 'edgeId' in pol2
+                  idIns1.push_back(std::pair<mcIdType,std::pair<bool,mcIdType> >(edgeId1,std::pair<bool,mcIdType>(descBg1[j]>0,offset1)));// it exists an edge into pol1 given by tuple (idIn1,direct1) that is colinear at edge 'edgeId' in pol2
                   //std::pair<edgeId1); direct1=descBg1[j]>0;
                 }
-              offset1+=intersectEdges1[edgeId1].size()/2;//offset1 is used to find the INTERP_KERNEL::Edge * instance into pol1 that will be part of edge into pol2
+              offset1+=ToIdType(intersectEdges1[edgeId1].size()/2);//offset1 is used to find the INTERP_KERNEL::Edge * instance into pol1 that will be part of edge into pol2
             }
           directos=idIns1.empty();
         }
@@ -534,25 +534,25 @@ void QuadraticPolygon::buildFromCrudeDataArray2(const std::map<int,INTERP_KERNEL
         }
       else
         {//there is subpart of edge 'edgeId' of pol2 inside pol1
-          const std::vector<int>& subEdge=intersectEdges2[edgeId];
+          const std::vector<mcIdType>& subEdge=intersectEdges2[edgeId];
           std::size_t nbOfSubEdges=subEdge.size()/2;
           for(std::size_t j=0;j<nbOfSubEdges;j++)
             {
-              int idBg=direct?subEdge[2*j]:subEdge[2*nbOfSubEdges-2*j-1];
-              int idEnd=direct?subEdge[2*j+1]:subEdge[2*nbOfSubEdges-2*j-2];
-              bool direction11,found=false;
-              bool direct1;//store if needed the direction in 1
-              int offset2;
-              std::size_t nbOfSubEdges1;
-              for(std::vector<std::pair<int,std::pair<bool,int> > >::const_iterator it=idIns1.begin();it!=idIns1.end() && !found;it++)
+              mcIdType idBg=direct?subEdge[2*j]:subEdge[2*nbOfSubEdges-2*j-1];
+              mcIdType idEnd=direct?subEdge[2*j+1]:subEdge[2*nbOfSubEdges-2*j-2];
+              bool direction11=false,found=false;
+              bool direct1=false;//store if needed the direction in 1
+              mcIdType offset2=0;
+              mcIdType nbOfSubEdges1=0;
+              for(std::vector<std::pair<mcIdType,std::pair<bool,mcIdType> > >::const_iterator it=idIns1.begin();it!=idIns1.end() && !found;it++)
                 {
-                  int idIn1=(*it).first;//store if needed the cell id in 1
+                  mcIdType idIn1=(*it).first;//store if needed the cell id in 1
                   direct1=(*it).second.first;
                   offset1=(*it).second.second;
-                  const std::vector<int>& subEdge1PossiblyAlreadyIn1=intersectEdges1[idIn1];
-                  nbOfSubEdges1=subEdge1PossiblyAlreadyIn1.size()/2;
+                  const std::vector<mcIdType>& subEdge1PossiblyAlreadyIn1=intersectEdges1[idIn1];
+                  nbOfSubEdges1=ToIdType(subEdge1PossiblyAlreadyIn1.size()/2);
                   offset2=0;
-                  for(std::size_t k=0;k<nbOfSubEdges1 && !found;k++)
+                  for(mcIdType k=0;k<nbOfSubEdges1 && !found;k++)
                     {//perform a loop on all subedges of pol1 that includes edge 'edgeId' of pol2. For the moment we iterate only on subedges of ['idIn1']... To improve
                       if(subEdge1PossiblyAlreadyIn1[2*k]==idBg && subEdge1PossiblyAlreadyIn1[2*k+1]==idEnd)
                         { direction11=true; found=true; }
@@ -573,7 +573,7 @@ void QuadraticPolygon::buildFromCrudeDataArray2(const std::map<int,INTERP_KERNEL
                 }
               else
                 {//the current subedge of edge 'edgeId' of pol2 is part of the colinear edge 'idIn1' of pol1 -> reuse Edge instance of pol1
-                  ElementaryEdge *e=pol1[offset1+(direct1?offset2:nbOfSubEdges1-offset2-1)];
+                  ElementaryEdge *e=pol1[FromIdType<int>(offset1+(direct1?offset2:nbOfSubEdges1-offset2-1))];
                   Edge *ee=e->getPtr();
                   ee->incrRef();
                   ElementaryEdge *e2=new ElementaryEdge(ee,!(direct1^direction11));
@@ -589,39 +589,39 @@ void QuadraticPolygon::buildFromCrudeDataArray2(const std::map<int,INTERP_KERNEL
  * Method expected to be called on pol2. Every params not suffixed by numbered are supposed to refer to pol2 (this).
  * Method to find edges that are ON.
  */
-void QuadraticPolygon::updateLocOfEdgeFromCrudeDataArray2(const int *descBg, const int *descEnd, const std::vector<std::vector<int> >& intersectEdges,
-                                                          const INTERP_KERNEL::QuadraticPolygon& pol1, const int *descBg1, const int *descEnd1,
-                                                          const std::vector<std::vector<int> >& intersectEdges1, const std::vector< std::vector<int> >& colinear1) const
+void QuadraticPolygon::updateLocOfEdgeFromCrudeDataArray2(const mcIdType *descBg, const mcIdType *descEnd, const std::vector<std::vector<mcIdType> >& intersectEdges,
+                                                          const INTERP_KERNEL::QuadraticPolygon& pol1, const mcIdType *descBg1, const mcIdType *descEnd1,
+                                                          const std::vector<std::vector<mcIdType> >& intersectEdges1, const std::vector< std::vector<mcIdType> >& colinear1) const
 {
   std::size_t nbOfSeg=std::distance(descBg,descEnd);
   for(std::size_t i=0;i<nbOfSeg;i++)//loop over all edges of pol2
     {
       bool direct=descBg[i]>0;
-      int edgeId=abs(descBg[i])-1;//current edge id of pol2
-      const std::vector<int>& c=colinear1[edgeId];
+      mcIdType edgeId=std::abs(descBg[i])-1;//current edge id of pol2
+      const std::vector<mcIdType>& c=colinear1[edgeId];
       if(c.empty())
         continue;
-      const std::vector<int>& subEdge=intersectEdges[edgeId];
+      const std::vector<mcIdType>& subEdge=intersectEdges[edgeId];
       std::size_t nbOfSubEdges=subEdge.size()/2;
       //
       std::size_t nbOfEdgesIn1=std::distance(descBg1,descEnd1);
-      int offset1=0;
+      mcIdType offset1=0;
       for(std::size_t j=0;j<nbOfEdgesIn1;j++)
         {
-          int edgeId1=abs(descBg1[j])-1;
+          mcIdType edgeId1=std::abs(descBg1[j])-1;
           if(std::find(c.begin(),c.end(),edgeId1)!=c.end())
             {
               for(std::size_t k=0;k<nbOfSubEdges;k++)
                 {
-                  int idBg=direct?subEdge[2*k]:subEdge[2*nbOfSubEdges-2*k-1];
-                  int idEnd=direct?subEdge[2*k+1]:subEdge[2*nbOfSubEdges-2*k-2];
-                  int idIn1=edgeId1;
+                  mcIdType idBg=direct?subEdge[2*k]:subEdge[2*nbOfSubEdges-2*k-1];
+                  mcIdType idEnd=direct?subEdge[2*k+1]:subEdge[2*nbOfSubEdges-2*k-2];
+                  mcIdType idIn1=edgeId1;
                   bool direct1=descBg1[j]>0;
-                  const std::vector<int>& subEdge1PossiblyAlreadyIn1=intersectEdges1[idIn1];
-                  std::size_t nbOfSubEdges1=subEdge1PossiblyAlreadyIn1.size()/2;
-                  int offset2=0;
+                  const std::vector<mcIdType>& subEdge1PossiblyAlreadyIn1=intersectEdges1[idIn1];
+                  mcIdType nbOfSubEdges1=ToIdType(subEdge1PossiblyAlreadyIn1.size()/2);
+                  mcIdType offset2=0;
                   bool found=false;
-                  for(std::size_t kk=0;kk<nbOfSubEdges1 && !found;kk++)
+                  for(mcIdType kk=0;kk<nbOfSubEdges1 && !found;kk++)
                     {
                       found=(subEdge1PossiblyAlreadyIn1[2*kk]==idBg && subEdge1PossiblyAlreadyIn1[2*kk+1]==idEnd) || (subEdge1PossiblyAlreadyIn1[2*kk]==idEnd && subEdge1PossiblyAlreadyIn1[2*kk+1]==idBg);
                       if(!found)
@@ -629,17 +629,17 @@ void QuadraticPolygon::updateLocOfEdgeFromCrudeDataArray2(const int *descBg, con
                     }
                   if(found)
                     {
-                      ElementaryEdge *e=pol1[offset1+(direct1?offset2:nbOfSubEdges1-offset2-1)];
+                      ElementaryEdge *e=pol1[(int)(offset1+(direct1?offset2:nbOfSubEdges1-offset2-1))];
                       e->getPtr()->declareOn();
                     }
                 }
             }
-          offset1+=intersectEdges1[edgeId1].size()/2;//offset1 is used to find the INTERP_KERNEL::Edge * instance into pol1 that will be part of edge into pol2
+          offset1+=ToIdType(intersectEdges1[edgeId1].size()/2);//offset1 is used to find the INTERP_KERNEL::Edge * instance into pol1 that will be part of edge into pol2
         }
     }
 }
 
-void QuadraticPolygon::appendCrudeData(const std::map<INTERP_KERNEL::Node *,int>& mapp, double xBary, double yBary, double fact, int offset, std::vector<double>& addCoordsQuadratic, std::vector<int>& conn, std::vector<int>& connI) const
+void QuadraticPolygon::appendCrudeData(const std::map<INTERP_KERNEL::Node *,mcIdType>& mapp, double xBary, double yBary, double fact, mcIdType offset, std::vector<double>& addCoordsQuadratic, std::vector<mcIdType>& conn, std::vector<mcIdType>& connI) const
 {
   int nbOfNodesInPg=0;
   bool presenceOfQuadratic=presenceOfQuadraticEdge();
@@ -648,14 +648,14 @@ void QuadraticPolygon::appendCrudeData(const std::map<INTERP_KERNEL::Node *,int>
     {
       Node *tmp=0;
       tmp=(*it)->getStartNode();
-      std::map<INTERP_KERNEL::Node *,int>::const_iterator it1=mapp.find(tmp);
+      std::map<INTERP_KERNEL::Node *,mcIdType>::const_iterator it1=mapp.find(tmp);
       conn.push_back((*it1).second);
       nbOfNodesInPg++;
     }
   if(presenceOfQuadratic)
     {
       int j=0;
-      int off=offset+((int)addCoordsQuadratic.size())/2;
+      mcIdType off=offset+ToIdType(addCoordsQuadratic.size())/2;
       for(std::list<ElementaryEdge *>::const_iterator it=_sub_edges.begin();it!=_sub_edges.end();it++,j++,nbOfNodesInPg++)
         {
           INTERP_KERNEL::Node *node=(*it)->getPtr()->buildRepresentantOfMySelf();
@@ -672,13 +672,13 @@ void QuadraticPolygon::appendCrudeData(const std::map<INTERP_KERNEL::Node *,int>
 /*!
  * This method make the hypothesis that \a this and \a other are split at the minimum into edges that are fully IN, OUT or ON.
  * This method returns newly created polygons in \a conn and \a connI and the corresponding ids ( \a idThis, \a idOther) are stored respectively into \a nbThis and \a nbOther.
- * @param [in,out] edgesThis, parameter that keep informed the caller about the edges in this not shared by the result of intersection of \a this with \a other
- * @param [in,out] edgesBoundaryOther, parameter that stores all edges in result of intersection that are not
+ * @param [in,out] edgesThis parameter that keep informed the caller about the edges in this not shared by the result of intersection of \a this with \a other
+ * @param [in,out] edgesBoundaryOther parameter that stores all edges in result of intersection that are not
  */
 void QuadraticPolygon::buildPartitionsAbs(QuadraticPolygon& other, std::set<INTERP_KERNEL::Edge *>& edgesThis, std::set<INTERP_KERNEL::Edge *>& edgesBoundaryOther,
-                                          const std::map<INTERP_KERNEL::Node *,int>& mapp, int idThis, int idOther, int offset,
-                                          std::vector<double>& addCoordsQuadratic, std::vector<int>& conn, std::vector<int>& connI,
-                                          std::vector<int>& nbThis, std::vector<int>& nbOther)
+                                          const std::map<INTERP_KERNEL::Node *,mcIdType>& mapp, mcIdType idThis, mcIdType idOther, mcIdType offset,
+                                          std::vector<double>& addCoordsQuadratic, std::vector<mcIdType>& conn, std::vector<mcIdType>& connI,
+                                          std::vector<mcIdType>& nbThis, std::vector<mcIdType>& nbOther)
 {
   double xBaryBB, yBaryBB;
   double fact=normalizeExt(&other, xBaryBB, yBaryBB);
@@ -1080,7 +1080,7 @@ std::list<QuadraticPolygon *> QuadraticPolygon::zipConsecutiveInSegments() const
 }
 
 /*!
- * @param [in] pol1zip is a list of set of edges (=an opened polygon) coming from split polygon 1.
+ * @param [in] pol1Zip is a list of set of edges (=an opened polygon) coming from split polygon 1.
  * @param [in] pol1 should be considered as pol1Simplified.
  * @param [in] pol2 is split pol2.
  * @param [out] results the resulting \b CLOSED polygons.
@@ -1089,7 +1089,8 @@ void QuadraticPolygon::ClosePolygons(std::list<QuadraticPolygon *>& pol1Zip, con
                                      std::vector<QuadraticPolygon *>& results)
 {
   bool directionKnownInPol2=false;
-  bool directionInPol2;
+  bool directionInPol2=false;
+  bool needCleaning = false;
   for(std::list<QuadraticPolygon *>::iterator iter=pol1Zip.begin();iter!=pol1Zip.end();)
     {
       // Build incrementally the full closed cells from the consecutive line parts already built in pol1Zip.
@@ -1097,14 +1098,17 @@ void QuadraticPolygon::ClosePolygons(std::list<QuadraticPolygon *>& pol1Zip, con
       // This process can produce several cells.
       if((*iter)->completed())
         {
+          if (needCleaning)
+            (*iter)->cleanDegeneratedConsecutiveEdges();
           results.push_back(*iter);
           directionKnownInPol2=false;
+          needCleaning=false;
           iter=pol1Zip.erase(iter);
           continue;
         }
       if(!directionKnownInPol2)
         {
-          if(!(*iter)->haveIAChanceToBeCompletedBy(pol1,pol2,directionInPol2))
+          if(!(*iter)->haveIAChanceToBeCompletedBy(pol1,pol2,directionInPol2, needCleaning))
             { delete *iter; iter=pol1Zip.erase(iter); continue; }
           else
             directionKnownInPol2=true;
@@ -1125,12 +1129,15 @@ void QuadraticPolygon::ClosePolygons(std::list<QuadraticPolygon *>& pol1Zip, con
 /*!
  * 'this' is expected to be set of edges (not closed) of pol1 split.
  */
-bool QuadraticPolygon::haveIAChanceToBeCompletedBy(const QuadraticPolygon& pol1NotSplitted, const QuadraticPolygon& pol2Splitted, bool& direction) const
+bool QuadraticPolygon::haveIAChanceToBeCompletedBy(const QuadraticPolygon& pol1NotSplitted, const QuadraticPolygon& pol2Splitted,
+                                                   bool& direction, bool& needCleaning) const
 {
+  needCleaning = false;
   IteratorOnComposedEdge it2(const_cast<QuadraticPolygon *>(&pol2Splitted));
   bool found=false;
   Node *n=getEndNode();
   ElementaryEdge *cur=it2.current();
+  // Find edge in pol2 whose start node is the end node of the current piece in pol1Zip (*this)
   for(it2.first();!it2.finished() && !found;)
     {
       cur=it2.current();
@@ -1146,20 +1153,26 @@ bool QuadraticPolygon::haveIAChanceToBeCompletedBy(const QuadraticPolygon& pol1N
     {
       if(e->getPtr()==cur->getPtr())
         {
-          direction=false;
-          it2.previousLoop();
+          // if we have the same edge, several possibilities:
+          //    - either this means that pol1 and pol2 have opposite orientation (since we matched end node with start node before)
+          //    - or (more tricky, see testIntersect2DMeshes11()) that pol1 and pol2 have same orientation but 'this' turns in such
+          //    a way that it attaches to pol2 on an edge in opposite orientation.
+          // To sort this out, inspect localisation of next edge in pol2 wrt pol1NotSplitted.
+          it2.nextLoop();
           cur=it2.current();
           Node *repr=cur->getPtr()->buildRepresentantOfMySelf();
           bool ret=pol1NotSplitted.isInOrOut(repr);
           repr->decrRef();
+          direction = ret;
+          needCleaning = ret; // if true we are in tricky case 2 above, we know that we will produce two consecutive overlapping edges in result
           return ret;
         }
-      else
+      else  // here we don't need to go prev or next:
         {
-          direction=true;
           Node *repr=cur->getPtr()->buildRepresentantOfMySelf();
           bool ret=pol1NotSplitted.isInOrOut(repr);
           repr->decrRef();
+          direction = ret;
           return ret;
         }
     }
@@ -1231,9 +1244,9 @@ std::list<QuadraticPolygon *>::iterator QuadraticPolygon::CheckInList(Node *n, s
 *      intersecting cells
 */
 void QuadraticPolygon::ComputeResidual(const QuadraticPolygon& pol1, const std::set<Edge *>& notUsedInPol1, const std::set<Edge *>& edgesInPol2OnBoundary,
-                                       const std::map<INTERP_KERNEL::Node *,int>& mapp, int offset, int idThis,
-                                       std::vector<double>& addCoordsQuadratic, std::vector<int>& conn,
-                                       std::vector<int>& connI, std::vector<int>& nb1, std::vector<int>& nb2)
+                                       const std::map<INTERP_KERNEL::Node *,mcIdType>& mapp, mcIdType offset, mcIdType idThis,
+                                       std::vector<double>& addCoordsQuadratic, std::vector<mcIdType>& conn,
+                                       std::vector<mcIdType>& connI, std::vector<mcIdType>& nb1, std::vector<mcIdType>& nb2)
 {
   // Initialise locations on pol1. Remember that edges found in 'notUsedInPol1' are also part of the edges forming pol1.
   pol1.initLocations();
@@ -1291,24 +1304,55 @@ void QuadraticPolygon::ComputeResidual(const QuadraticPolygon& pol1, const std::
       // retPolsUnderConstruction initially empty -> see if(!pol1Zip.empty()) below ...
       for(std::list<QuadraticPolygon *>::iterator itConstr=retPolsUnderContruction.begin();itConstr!=retPolsUnderContruction.end();)
         {
-          if((*itConstr)->getStartNode()==(*itConstr)->getEndNode())  // reconstruction of a cell is finished
-            {
-              itConstr++;
-              continue;
-            }
-          Node *curN=(*itConstr)->getEndNode();
-          bool smthHappened=false;
-          // Complete a partially reconstructed polygon with boundary edges by matching nodes:
+          Node *startN = (*itConstr)->getStartNode();
+          Node *curN = (*itConstr)->getEndNode();
+          if(startN == curN)  // reconstruction of a cell is finished
+            {  itConstr++; continue;  }
+
+          bool smthHappened=false, doneEarly=false;
+          // Complete a partially reconstructed polygon with boundary edges of pol2 by matching nodes:
           for(std::list<Edge *>::iterator it2=edgesInPol2OnBoundaryL.begin();it2!=edgesInPol2OnBoundaryL.end();)
             {
-              if(curN==(*it2)->getStartNode())
-                { (*it2)->incrRef(); (*itConstr)->pushBack(new ElementaryEdge(*it2,true)); curN=(*it2)->getEndNode(); smthHappened=true; it2=edgesInPol2OnBoundaryL.erase(it2); }
-              else if(curN==(*it2)->getEndNode())
-                { (*it2)->incrRef(); (*itConstr)->pushBack(new ElementaryEdge(*it2,false)); curN=(*it2)->getStartNode(); smthHappened=true; it2=edgesInPol2OnBoundaryL.erase(it2); }
+              if(curN==(*it2)->getEndNode())  // only end node should be considered if orientation is correct for input meshes
+                                              // in the funny case of cells exactly included (see test case testIntersect2DMeshesTmp13) this is mandatory to take edges from pol2 in the right order.
+                {
+                  (*it2)->incrRef();
+                  (*itConstr)->pushBack(new ElementaryEdge(*it2,false));
+                  curN=(*it2)->getStartNode();
+                  smthHappened=true;
+                  it2=edgesInPol2OnBoundaryL.erase(it2);
+                }
               else
                 it2++;
+              if (curN == startN) // we might end here
+                { doneEarly = true;  break;  }
             }
-          if(smthHappened)
+
+          // It might be the case that the lookup on start nodes made above failed because pol2 is wrongly oriented.
+          // Be somewhat flexible and keep on supporting this case here (useful for voronisation notably):
+          if(!smthHappened)
+            {
+              for(std::list<Edge *>::iterator it2=edgesInPol2OnBoundaryL.begin();it2!=edgesInPol2OnBoundaryL.end();)
+                {
+                  if(curN==(*it2)->getStartNode())
+                    {
+                      (*it2)->incrRef();
+                      (*itConstr)->pushBack(new ElementaryEdge(*it2,true));
+                      curN=(*it2)->getEndNode();
+                      smthHappened=true;
+                      it2=edgesInPol2OnBoundaryL.erase(it2);
+                    }
+                  else
+                    it2++;
+                  if (curN == startN) // we might end here
+                      { doneEarly = true;  break;  }
+                }
+            }
+
+          if (doneEarly)
+            {  itConstr++; continue;  }
+
+          if(smthHappened) // Now continue the construction by finding the next bit in pol1Zip. Not too sure what are the cases where the boolean if False ...
             {
               for(std::list<QuadraticPolygon *>::iterator itZip=pol1Zip.begin();itZip!=pol1Zip.end();)
                 {
@@ -1325,7 +1369,7 @@ void QuadraticPolygon::ComputeResidual(const QuadraticPolygon& pol1, const std::
                     itZip++;
                 }
             }
-          else
+          else // Nothing happened.
             {
               for(std::list<ElementaryEdge *>::const_iterator it5=(*itConstr)->_sub_edges.begin();it5!=(*itConstr)->_sub_edges.end();it5++)
                 {
@@ -1363,8 +1407,9 @@ void QuadraticPolygon::ComputeResidual(const QuadraticPolygon& pol1, const std::
   // Convert to integer connectivity:
   for(std::list<QuadraticPolygon *>::iterator itConstr=retPolsUnderContruction.begin();itConstr!=retPolsUnderContruction.end();itConstr++)
     {
-      if((*itConstr)->getStartNode()==(*itConstr)->getEndNode())  // take only fully closed reconstructed polygon (?? might there be others??)
+      if((*itConstr)->getStartNode()==(*itConstr)->getEndNode())  // take only fully closed reconstructed polygon
         {
+          (*itConstr)->cleanDegeneratedConsecutiveEdges();
           (*itConstr)->appendCrudeData(mapp,0.,0.,1.,offset,addCoordsQuadratic,conn,connI); nb1.push_back(idThis); nb2.push_back(-1);
           for(std::list<QuadraticPolygon *>::iterator it6=pol1ZipConsumed[*itConstr].begin();it6!=pol1ZipConsumed[*itConstr].end();it6++)
             delete *it6;