From 74d10b312f6ea6f61f5f5fbc61c0af07538a9ce6 Mon Sep 17 00:00:00 2001 From: abn Date: Wed, 6 Mar 2019 15:00:44 +0100 Subject: [PATCH] Intersect2DMeshes bug fix: mishandling of orientation in reconstruction. --- .../InterpKernelGeo2DQuadraticPolygon.cxx | 25 +++++++++++++----- .../InterpKernelGeo2DQuadraticPolygon.hxx | 3 ++- .../MEDCouplingIntersectTest.py | 26 +++++++++++++++++++ 3 files changed, 47 insertions(+), 7 deletions(-) diff --git a/src/INTERP_KERNEL/Geometric2D/InterpKernelGeo2DQuadraticPolygon.cxx b/src/INTERP_KERNEL/Geometric2D/InterpKernelGeo2DQuadraticPolygon.cxx index e10c446b3..3e4c6c10d 100644 --- a/src/INTERP_KERNEL/Geometric2D/InterpKernelGeo2DQuadraticPolygon.cxx +++ b/src/INTERP_KERNEL/Geometric2D/InterpKernelGeo2DQuadraticPolygon.cxx @@ -1090,6 +1090,7 @@ void QuadraticPolygon::ClosePolygons(std::list& pol1Zip, con { bool directionKnownInPol2=false; bool directionInPol2; + bool needCleaning = false; for(std::list::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& 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& 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(&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; } } diff --git a/src/INTERP_KERNEL/Geometric2D/InterpKernelGeo2DQuadraticPolygon.hxx b/src/INTERP_KERNEL/Geometric2D/InterpKernelGeo2DQuadraticPolygon.hxx index 7b063aafc..4814be460 100644 --- a/src/INTERP_KERNEL/Geometric2D/InterpKernelGeo2DQuadraticPolygon.hxx +++ b/src/INTERP_KERNEL/Geometric2D/InterpKernelGeo2DQuadraticPolygon.hxx @@ -93,7 +93,8 @@ namespace INTERP_KERNEL INTERPKERNEL_EXPORT void performLocatingOperationSlow(QuadraticPolygon& pol2) const; INTERPKERNEL_EXPORT static void SplitPolygonsEachOther(QuadraticPolygon& pol1, QuadraticPolygon& pol2, int& nbOfSplits); INTERPKERNEL_EXPORT std::vector buildIntersectionPolygons(const QuadraticPolygon& pol1, const QuadraticPolygon& pol2) const; - INTERPKERNEL_EXPORT bool haveIAChanceToBeCompletedBy(const QuadraticPolygon& pol1NotSplitted, const QuadraticPolygon& pol2Splitted, bool& direction) const; + INTERPKERNEL_EXPORT bool haveIAChanceToBeCompletedBy(const QuadraticPolygon& pol1NotSplitted, const QuadraticPolygon& pol2Splitted, + bool& direction, bool& needCleaning) const; INTERPKERNEL_EXPORT static void ComputeResidual(const QuadraticPolygon& pol1, const std::set& notUsedInPol1, const std::set& edgesInPol2OnBoundary, const std::map& mapp, int offset, int idThis, std::vector& addCoordsQuadratic, std::vector& conn, std::vector& connI, std::vector& nb1, std::vector& nb2); INTERPKERNEL_EXPORT void cleanDegeneratedConsecutiveEdges(); diff --git a/src/MEDCoupling_Swig/MEDCouplingIntersectTest.py b/src/MEDCoupling_Swig/MEDCouplingIntersectTest.py index e02005e82..bc9a30b00 100644 --- a/src/MEDCoupling_Swig/MEDCouplingIntersectTest.py +++ b/src/MEDCoupling_Swig/MEDCouplingIntersectTest.py @@ -439,6 +439,32 @@ class MEDCouplingIntersectTest(unittest.TestCase): self.assertEqual(e2, res2Tool.getValues()) pass + def testIntersect2DMeshes11(self): + """ Dealing properly with respective polygon orientation in QuadraticPolygon::haveIAChanceToBeCompletedBy() + The two polygons below have same orientation, but one edge of pol1 is colinear to pol2 in opposite directions. + """ + eps = 1.0e-8 + back = MEDCouplingUMesh('lback', 2) + coo = DataArrayDouble([(-2.5,-2.5),(-2.5,2.5),(2.5,2.5),(2.5,-2.5),(0,0),(0,1.66667),(1.66667,1.66667),(1.66667,0)]) + back.setCoords(coo) + c = DataArrayInt([5, 6, 7, 4, 5, 2, 3, 7, 6, 5, 3, 0, 1, 2, 6, 4, 7]) + cI = DataArrayInt([0, 4, 9, 17]) + back.setConnectivity(c, cI) + + tool = MEDCouplingUMesh('ltool', 2) + coo = DataArrayDouble([(0,0),(0,2.5),(2.5,2.5),(2.5,0)]) + tool.setCoords(coo) + c = DataArrayInt([5, 0, 1, 2, 3]) + cI = DataArrayInt([0, 5]) + tool.setConnectivity(c, cI) + + result, res2Back, res2Tool = MEDCouplingUMesh.Intersect2DMeshes(back, tool, eps) + self.assertEqual(result.getNodalConnectivity().getValues(), [5, 7, 8, 6, 5, 7, 6, 10, 11, 5, 11, 3, 7, 5, 9, 10, 6, 8, 5, 8, 7, 3, 0, 1, 9]) + self.assertEqual(result.getNodalConnectivityIndex().getValues(), [0, 4, 9, 13, 18, 25]) + self.assertEqual(res2Back.getValues(), [0, 1, 1, 2, 2]) + self.assertEqual(res2Tool.getValues(), [0, 0, -1, 0, -1]) + pass + def testSwig2Intersect2DMeshesQuadra1(self): import cmath def createDiagCircle(lX, lY, R, cells=[0,1]): -- 2.39.2