From: abn Date: Mon, 6 Aug 2018 08:39:51 +0000 (+0200) Subject: [backport] New functionality: colinearizeKeepingConform2D(). X-Git-Url: http://git.salome-platform.org/gitweb/?a=commitdiff_plain;h=4a67ed435810cb1c88ba946d76649b18d2218bb3;p=tools%2Fmedcoupling.git [backport] New functionality: colinearizeKeepingConform2D(). --- diff --git a/src/INTERP_KERNEL/CellModel.cxx b/src/INTERP_KERNEL/CellModel.cxx index bd676e424..974e268df 100644 --- a/src/INTERP_KERNEL/CellModel.cxx +++ b/src/INTERP_KERNEL/CellModel.cxx @@ -529,7 +529,7 @@ namespace INTERP_KERNEL sonNodalConn[1]=nodalConn[(sonId+1)%lgth]; return 2; } - else + else // NORM_QPOLYG { sonNodalConn[0]=nodalConn[sonId]; sonNodalConn[1]=nodalConn[(sonId+1)%(lgth/2)]; diff --git a/src/MEDCoupling/MEDCouplingUMesh.cxx b/src/MEDCoupling/MEDCouplingUMesh.cxx index 6db248d97..31a0d6e6d 100644 --- a/src/MEDCoupling/MEDCouplingUMesh.cxx +++ b/src/MEDCoupling/MEDCouplingUMesh.cxx @@ -10518,6 +10518,27 @@ DataArrayInt *MEDCouplingUMesh::conformize2D(double eps) * \sa MEDCouplingUMesh::conformize2D, MEDCouplingUMesh::mergeNodes, MEDCouplingUMesh::convexEnvelop2D. */ DataArrayInt *MEDCouplingUMesh::colinearize2D(double eps) +{ + internalColinearize2D(eps, false); +} + +/*! + * Performs exactly the same job as colinearize2D, except that this function does not create new non-conformal cells. + * In a given 2D cell, if two edges are colinear and the junction point is used by a third edge, the two edges will not be + * merged, contrary to colinearize2D(). + * + * \sa MEDCouplingUMesh::colinearize2D + */ +DataArrayInt *MEDCouplingUMesh::colinearizeKeepingConform2D(double eps) +{ + internalColinearize2D(eps, true); +} + + +/*! + * \param stayConform is set to True, will not fuse two edges sharing a node that has (strictly) more than 2 egdes connected to it + */ +DataArrayInt *MEDCouplingUMesh::internalColinearize2D(double eps, bool stayConform) { MEDCouplingAutoRefCountObjectPtr ret(DataArrayInt::New()); ret->alloc(0,1); checkCoherency(); @@ -10528,12 +10549,27 @@ DataArrayInt *MEDCouplingUMesh::colinearize2D(double eps) int nbOfCells(getNumberOfCells()),nbOfNodes(getNumberOfNodes()); const int *cptr(_nodal_connec->begin()),*ciptr(_nodal_connec_index->begin()); MEDCouplingAutoRefCountObjectPtr newc(DataArrayInt::New()),newci(DataArrayInt::New()); newci->alloc(nbOfCells+1,1); newc->alloc(0,1); newci->setIJ(0,0,0); + std::map forbiddenPoints; // list of points that can not be removed (or it will break conformity) + if(stayConform) // + { + // A point that is used by more than 2 edges can not be removed without breaking conformity: + MEDCouplingAutoRefCountObjectPtr desc1(DataArrayInt::New()),descI1(DataArrayInt::New()),revDesc1(DataArrayInt::New()),revDescI1(DataArrayInt::New()); + MEDCouplingAutoRefCountObjectPtr mDesc1D(buildDescendingConnectivity(desc1,descI1,revDesc1,revDescI1)); + MEDCouplingAutoRefCountObjectPtr desc2(DataArrayInt::New()),descI2(DataArrayInt::New()),revDesc2(DataArrayInt::New()),revDescI2(DataArrayInt::New()); + MEDCouplingAutoRefCountObjectPtr mDesc0D(mDesc1D->buildDescendingConnectivity(desc2,descI2,revDesc2,revDescI2)); + MEDCouplingAutoRefCountObjectPtr dsi(revDescI2->deltaShiftIndex()); + MEDCouplingAutoRefCountObjectPtr ids(dsi->getIdsInRange(3, 10000)); + const int * cPtr(mDesc0D->getNodalConnectivity()->begin()); + for(const int* it = ids->begin(); it != ids->end(); it++) + forbiddenPoints[cPtr[2*(*it)+1]] = true; // we know that a 0D mesh has a connectivity of the form [NORM_POINT1, i1, NORM_POINT1, i2, ...] + } + MEDCouplingAutoRefCountObjectPtr appendedCoords(DataArrayDouble::New()); appendedCoords->alloc(0,1);//1 not 2 it is not a bug. const double *coords(_coords->begin()); int *newciptr(newci->getPointer()); for(int i=0;ipushBackSilent(i); newciptr[1]=newc->getNumberOfTuples(); } @@ -11768,13 +11804,16 @@ void EnterTheResultOf2DCellEnd(const INTERP_KERNEL::Edge *e, int start, int stp, * Returns true if a colinearization has been found in the given cell. If false is returned the content pushed in \a newConnOfCell is equal to [ \a connBg , \a connEnd ) . * \a appendedCoords is a DataArrayDouble instance with number of components equal to one (even if the items are pushed by pair). */ -bool MEDCouplingUMesh::Colinearize2DCell(const double *coords, const int *connBg, const int *connEnd, int offset, DataArrayInt *newConnOfCell, DataArrayDouble *appendedCoords) +bool MEDCouplingUMesh::Colinearize2DCell(const double *coords, const int *connBg, const int *connEnd, int offset, + const std::map& forbiddenPoints, + DataArrayInt *newConnOfCell, DataArrayDouble *appendedCoords) { std::size_t sz(std::distance(connBg,connEnd)); if(sz<3)//3 because 2+1(for the cell type) and 2 is the minimal number of edges of 2D cell. throw INTERP_KERNEL::Exception("MEDCouplingUMesh::Colinearize2DCell : the input cell has invalid format !"); sz--; INTERP_KERNEL::AutoPtr tmpConn(new int[sz]); + INTERP_KERNEL::AutoPtr tmpConn2(new int[sz]); const INTERP_KERNEL::CellModel& cm(INTERP_KERNEL::CellModel::GetCellModel((INTERP_KERNEL::NormalizedCellType)connBg[0])); unsigned nbs(cm.getNumberOfSons2(connBg+1,sz)); unsigned nbOfHit(0); // number of fusions operated @@ -11794,10 +11833,15 @@ bool MEDCouplingUMesh::Colinearize2DCell(const double *coords, const int *connBg // This initializes posBaseElt. if(nbOfTurn==0) { - for(unsigned i=1;i::const_iterator itE(forbiddenPoints.end()); + if (forbiddenPoints.find(commPoint) != itE) // is the junction point in the list of points we can not remove? + break; + INTERP_KERNEL::Edge *eCand(MEDCouplingUMeshBuildQPFromEdge2(typeOfSon,tmpConn2,coords,m)); INTERP_KERNEL::EdgeIntersector *eint(INTERP_KERNEL::Edge::BuildIntersectorWith(e,eCand)); bool isColinear=eint->areColinears(); if(isColinear) @@ -11810,14 +11854,21 @@ bool MEDCouplingUMesh::Colinearize2DCell(const double *coords, const int *connBg eCand->decrRef(); if(!isColinear) break; + // Update last connectivity + std::copy((int *)tmpConn2, tmpConn2+sz, (int *)tmpConn); } } // Now move forward: const unsigned fwdStart = (nbOfTurn == 0 ? 0 : posBaseElt); // the first element to be inspected going forward - for(unsigned j=fwdStart+1;j::const_iterator itE(forbiddenPoints.end()); + if (forbiddenPoints.find(commPoint) != itE) // is the junction point in the list of points we can not remove? + break; + INTERP_KERNEL::Edge *eCand(MEDCouplingUMeshBuildQPFromEdge2(typeOfSon,tmpConn2,coords,m)); INTERP_KERNEL::EdgeIntersector *eint(INTERP_KERNEL::Edge::BuildIntersectorWith(e,eCand)); bool isColinear(eint->areColinears()); if(isColinear) @@ -11830,11 +11881,13 @@ bool MEDCouplingUMesh::Colinearize2DCell(const double *coords, const int *connBg eCand->decrRef(); if(!isColinear) break; + // Update last connectivity + std::copy((int *)tmpConn2, tmpConn2+sz, (int *)tmpConn); } //push [posBaseElt,posEndElt) in newConnOfCell using e - // The if clauses below are (volontary) not mutually exclusive: on a quad cell with 2 edges, the end of the connectivity is also its begining! + // The if clauses below are (voluntary) not mutually exclusive: on a quad cell with 2 edges, the end of the connectivity is also its beginning! if(nbOfTurn==0) - // at the begining of the connectivity (insert type) + // at the beginning of the connectivity (insert type) EnterTheResultOf2DCellFirst(e,posBaseElt,posEndElt,(int)nbs,cm.isQuadratic(),coords,connBg+1,offset,newConnOfCell,appendedCoords,middles); else if((nbOfHit+nbOfTurn) != (nbs-1)) // in the middle diff --git a/src/MEDCoupling/MEDCouplingUMesh.hxx b/src/MEDCoupling/MEDCouplingUMesh.hxx index 007c3b4ba..ab2a711f1 100644 --- a/src/MEDCoupling/MEDCouplingUMesh.hxx +++ b/src/MEDCoupling/MEDCouplingUMesh.hxx @@ -229,6 +229,7 @@ namespace ParaMEDMEM MEDCOUPLING_EXPORT DataArrayDouble *computePlaneEquationOf3DFaces() const; MEDCOUPLING_EXPORT DataArrayInt *conformize2D(double eps); MEDCOUPLING_EXPORT DataArrayInt *colinearize2D(double eps); + MEDCOUPLING_EXPORT DataArrayInt *colinearizeKeepingConform2D(double eps); MEDCOUPLING_EXPORT int split2DCells(const DataArrayInt *desc, const DataArrayInt *descI, const DataArrayInt *subNodesInSeg, const DataArrayInt *subNodesInSegI, const DataArrayInt *midOpt=0, const DataArrayInt *midOptI=0); MEDCOUPLING_EXPORT static MEDCouplingUMesh *Build0DMeshFromCoords(DataArrayDouble *da); MEDCOUPLING_EXPORT static MEDCouplingUMesh *MergeUMeshes(const MEDCouplingUMesh *mesh1, const MEDCouplingUMesh *mesh2); @@ -338,8 +339,9 @@ namespace ParaMEDMEM const int *desc, const int *descIndx, DataArrayInt *nodalRes, DataArrayInt *nodalResIndx, DataArrayInt *cellIds) const; void split2DCellsLinear(const DataArrayInt *desc, const DataArrayInt *descI, const DataArrayInt *subNodesInSeg, const DataArrayInt *subNodesInSegI); int split2DCellsQuadratic(const DataArrayInt *desc, const DataArrayInt *descI, const DataArrayInt *subNodesInSeg, const DataArrayInt *subNodesInSegI, const DataArrayInt *mid, const DataArrayInt *midI); - static bool Colinearize2DCell(const double *coords, const int *connBg, const int *connEnd, int offset, DataArrayInt *newConnOfCell, DataArrayDouble *appendedCoords); + static bool Colinearize2DCell(const double *coords, const int *connBg, const int *connEnd, int offset, const std::map& forbiddenPoints, DataArrayInt *newConnOfCell, DataArrayDouble *appendedCoords); static void ComputeAllTypesInternal(std::set& types, const DataArrayInt *nodalConnec, const DataArrayInt *nodalConnecIndex); + DataArrayInt *internalColinearize2D(double eps, bool stayConform); public: MEDCOUPLING_EXPORT static DataArrayInt *ComputeRangesFromTypeDistribution(const std::vector& code); MEDCOUPLING_EXPORT static const int N_MEDMEM_ORDER=24; diff --git a/src/MEDCoupling_Swig/MEDCouplingBasicsTest5.py b/src/MEDCoupling_Swig/MEDCouplingBasicsTest5.py index b9ba43ff6..f5611c09e 100644 --- a/src/MEDCoupling_Swig/MEDCouplingBasicsTest5.py +++ b/src/MEDCoupling_Swig/MEDCouplingBasicsTest5.py @@ -2118,6 +2118,26 @@ class MEDCouplingBasicsTest5(unittest.TestCase): self.assertTrue(m.getNodalConnectivityIndex().isEqual(DataArrayInt([0,7]))) pass + def testSwig2ColinearizeKeepingConform2D1(self): + eps = 1.0e-6 + # Just to get a nice coords array ... + mm = MEDCouplingCMesh(); arr = DataArrayDouble([0.0, 1.0,2.0]) + mm.setCoords(arr, arr); mm = mm.buildUnstructured(); coo = mm.getCoords() + + mesh = MEDCouplingUMesh("M", 2) + mesh.setCoords(coo) + c = [NORM_POLYGON, 0,1,4,7,6,3, NORM_QUAD4, 1,2,5,4, NORM_QUAD4,4,5,8,7] + cI = [0, 7,12,17] + mm.setConnectivity(DataArrayInt(c),DataArrayInt(cI)) + mm.checkConsistencyLight() + + mm.colinearizeKeepingConform2D(eps) + c = mm.getNodalConnectivity().getValues() + cI = mm.getNodalConnectivityIndex().getValues() + self.assertEqual(c, [NORM_POLYGON, 0, 1, 4, 7, 6, NORM_POLYGON, 1, 2, 5, 4, NORM_POLYGON, 4, 5, 8, 7]) + self.assertEqual(cI, [0,6,11,16]) + pass + def testSwig2BoundingBoxForBBTree1(self): """ This test appears simple but it checks that bounding box are correctly computed for quadratic polygons. It can help a lot to reduce the amount of intersections ! """ diff --git a/src/MEDCoupling_Swig/MEDCouplingCommon.i b/src/MEDCoupling_Swig/MEDCouplingCommon.i index 582d22ada..23ef1f468 100644 --- a/src/MEDCoupling_Swig/MEDCouplingCommon.i +++ b/src/MEDCoupling_Swig/MEDCouplingCommon.i @@ -285,6 +285,7 @@ using namespace INTERP_KERNEL; %newobject ParaMEDMEM::MEDCouplingUMesh::buildNewNumberingFromCommNodesFrmt; %newobject ParaMEDMEM::MEDCouplingUMesh::conformize2D; %newobject ParaMEDMEM::MEDCouplingUMesh::colinearize2D; +%newobject ParaMEDMEM::MEDCouplingUMesh::colinearizeKeepingConform2D; %newobject ParaMEDMEM::MEDCouplingUMesh::rearrange2ConsecutiveCellTypes; %newobject ParaMEDMEM::MEDCouplingUMesh::sortCellsInMEDFileFrmt; %newobject ParaMEDMEM::MEDCouplingUMesh::getRenumArrForMEDFileFrmt; @@ -1773,6 +1774,7 @@ namespace ParaMEDMEM //tools DataArrayInt *conformize2D(double eps) throw(INTERP_KERNEL::Exception); DataArrayInt *colinearize2D(double eps) throw(INTERP_KERNEL::Exception); + DataArrayInt *colinearizeKeepingConform2D(double eps) throw(INTERP_KERNEL::Exception); void shiftNodeNumbersInConn(int delta) throw(INTERP_KERNEL::Exception); std::vector getQuadraticStatus() const throw(INTERP_KERNEL::Exception); DataArrayInt *findCellIdsOnBoundary() const throw(INTERP_KERNEL::Exception);