From c5de20b21dd06a88663b9c6b8729f6749a61a812 Mon Sep 17 00:00:00 2001 From: Anida Khizar Date: Thu, 17 Feb 2022 22:23:54 +0100 Subject: [PATCH] [Bug fix] : creation of redundant nodes in BuildSlice3D --- src/MEDCoupling/MEDCouplingUMesh.cxx | 6 +++- src/MEDCoupling/MEDCouplingUMesh_internal.cxx | 15 +++++++-- .../MEDCouplingBasicsTest3.py | 33 ++++++++++++++++++- 3 files changed, 49 insertions(+), 5 deletions(-) diff --git a/src/MEDCoupling/MEDCouplingUMesh.cxx b/src/MEDCoupling/MEDCouplingUMesh.cxx index a9181f3e7..0bbeeef39 100755 --- a/src/MEDCoupling/MEDCouplingUMesh.cxx +++ b/src/MEDCoupling/MEDCouplingUMesh.cxx @@ -3838,7 +3838,11 @@ MEDCouplingUMesh *MEDCouplingUMesh::buildSlice3D(const double *origin, const dou revDesc2=0; revDescIndx2=0; MCAuto mDesc1=mDesc2->buildDescendingConnectivity(desc1,descIndx1,revDesc1,revDescIndx1);//meshDim==1 spaceDim==3 revDesc1=0; revDescIndx1=0; - mDesc1->fillCellIdsToKeepFromNodeIds(&nodes[0],&nodes[0]+nodes.size(),true,cellIds1D); + //Marking all 1D cells that contained at least one node located on the plane + //the intersection between those cells and the plane, which consist of the nodes previously tagged, thus don't need to be computed afterwards + //(if said intersection is computed in MEDCouplingUMesh::split3DCurveWithPlane, then we might create additional nodes + //due to accuracy errors when the needed nodes already exist) + mDesc1->fillCellIdsToKeepFromNodeIds(&nodes[0],&nodes[0]+nodes.size(),false,cellIds1D); MCAuto cellIds1DTmp(cellIds1D); // std::vector cut3DCurve(mDesc1->getNumberOfCells(),-2); diff --git a/src/MEDCoupling/MEDCouplingUMesh_internal.cxx b/src/MEDCoupling/MEDCouplingUMesh_internal.cxx index 7bbe37f8e..7708548e0 100755 --- a/src/MEDCoupling/MEDCouplingUMesh_internal.cxx +++ b/src/MEDCoupling/MEDCouplingUMesh_internal.cxx @@ -547,7 +547,7 @@ void MEDCouplingUMesh::fillCellIdsToKeepFromNodeIds(const mcIdType *begin, const * This method has one in/out parameter : 'cut3DCurve'. * Param 'cut3DCurve' is expected to be of size 'this->getNumberOfCells()'. For each i in [0,'this->getNumberOfCells()') * if cut3DCurve[i]==-2, it means that for cell #i in \a this nothing has been detected previously. - * if cut3DCurve[i]==-1, it means that cell#i has been already detected to be fully part of plane defined by ('origin','vec'). + * if cut3DCurve[i]==-1, it means that cell#i has been already detected to be fully (or partially) 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& cut3DCurve) @@ -1624,11 +1624,20 @@ void MEDCouplingUMesh::AssemblyForSplitFrom3DCurve(const std::vector& res.push_back(status); else { - res.push_back(nodal3DCurve[nodalIndx3DCurve[edgeId]+1]); - res.push_back(nodal3DCurve[nodalIndx3DCurve[edgeId]+2]); + const mcIdType& node1 = nodal3DCurve[nodalIndx3DCurve[edgeId]+1]; + const mcIdType& node2 = nodal3DCurve[nodalIndx3DCurve[edgeId]+2]; + // Here, we have an edge that has either one or both of its nodes intersecting the plane + // we're only adding the nodes from the edges fully contained in the plane + if(std::find(nodesOnP.begin(), nodesOnP.end(),node1) != nodesOnP.end() + && std::find(nodesOnP.begin(), nodesOnP.end(),node2) != nodesOnP.end()) + { + res.push_back(node1); + res.push_back(node2); + } } } } + switch(res.size()) { case 2: diff --git a/src/MEDCoupling_Swig/MEDCouplingBasicsTest3.py b/src/MEDCoupling_Swig/MEDCouplingBasicsTest3.py index 521d00736..5c56013b9 100644 --- a/src/MEDCoupling_Swig/MEDCouplingBasicsTest3.py +++ b/src/MEDCoupling_Swig/MEDCouplingBasicsTest3.py @@ -2608,8 +2608,39 @@ class MEDCouplingBasicsTest3(unittest.TestCase): for i in range(135): self.assertAlmostEqual(expected11[i],slice1.getCoords().getIJ(0,i),12); pass + ## + coo = DataArrayDouble([1.8686305176182501,0.85678785370447097,7.4, 1.8686305176182501,0.85678785370447097,7.5, 1.8663408654534299,0.82910039403216995,7.4, 1.8663408654534299,0.82910039403216995,7.5, 1.8370211426271501,0.83286926189135702,7.4, 1.8370211426271501,0.83286926189135702,7.5, 1.84595839792064,0.86012397150595199,7.4, 1.84595839792064,0.86012397150595199,7.5], 8,3) + conn = DataArrayInt([18,0,2,4,6,1,3,5,7]) + connI = DataArrayInt([0,9]) + + mesh3D_2=MEDCouplingUMesh.New(); + mesh3D_2.setName("3DMesh_2"); + mesh3D_2.setMeshDimension(3); + mesh3D_2.setCoords(coo); + mesh3D_2.setConnectivity(conn,connI,True); + + expected12=[0] + expected13=[5,5,9,8,4] + expected14=[0,5] + expected15=[1.8686305176182501,0.85678785370447097,7.4, 1.8686305176182501,0.85678785370447097,7.5, 1.8663408654534299,0.82910039403216995,7.4, 1.8663408654534299,0.82910039403216995,7.5, 1.8370211426271501,0.83286926189135702,7.4, 1.8370211426271501,0.83286926189135702,7.5, 1.84595839792064,0.86012397150595199, 7.4, 1.84595839792064,0.86012397150595199,7.5, 1.8666525378798646,0.83286927118760312,7.4, 1.8666525378798646,0.83286927118760312,7.5] + + y_cut = 0.8328692711876031 + ori, vec = [0.0, y_cut, 0.0], [0.0,1.0,0.0] + slice1, ids = mesh3D_2.buildSlice3D(DataArrayDouble(ori,1,3), DataArrayDouble(vec,1,3), 1.0e-8) + self.assertEqual(2,slice1.getMeshDimension()); + self.assertEqual(3,slice1.getSpaceDimension()); + self.assertEqual(10,slice1.getNumberOfNodes()); + self.assertEqual(1,slice1.getNumberOfCells()); + self.assertEqual(1,ids.getNumberOfTuples()); + self.assertEqual(5,slice1.getNodalConnectivity().getNumberOfTuples()); + self.assertEqual(2,slice1.getNodalConnectivityIndex().getNumberOfTuples()); + self.assertEqual(expected12,ids.getValues()); + self.assertEqual(expected13,slice1.getNodalConnectivity().getValues()); + self.assertEqual(expected14,slice1.getNodalConnectivityIndex().getValues()); + for i in range(27): + self.assertAlmostEqual(expected15[i],slice1.getCoords().getIJ(0,i),12); + pass pass - def testBuildSlice3DSurf1(self): mesh3D,mesh2D=MEDCouplingDataForTest.build3DExtrudedUMesh_1(); mesh2D=mesh3D.buildDescendingConnectivity()[0]; -- 2.39.2