From 6a04d0269feb495fc286cd0fa7c39ddf91f931ef Mon Sep 17 00:00:00 2001 From: eap Date: Mon, 8 Nov 2021 19:28:43 +0300 Subject: [PATCH] bos #26452 [EDF] (2021) SMESH: orientation of faces --- .../examples/medcouplingexamplesmeshes.doxy | 5 +- .../doxygen/fakesources/MEDCouplingUMesh.C | 1 + doc/user/input/data_analysis.rst | 12 + src/MEDCoupling/MEDCouplingUMesh.cxx | 408 ++++++++++++++++++ src/MEDCoupling/MEDCouplingUMesh.hxx | 1 + .../Test/MEDCouplingExamplesTest.cxx | 11 + .../MEDCouplingBasicsTest1.py | 63 +++ src/MEDCoupling_Swig/MEDCouplingCommon.i | 1 + .../MEDCouplingExamplesTest.py | 11 + .../UsersGuideExamplesTest.py | 13 + 10 files changed, 525 insertions(+), 1 deletion(-) diff --git a/doc/developer/doxygen/doxfiles/examples/medcouplingexamplesmeshes.doxy b/doc/developer/doxygen/doxfiles/examples/medcouplingexamplesmeshes.doxy index ffa1d4112..491e4f01d 100644 --- a/doc/developer/doxygen/doxfiles/examples/medcouplingexamplesmeshes.doxy +++ b/doc/developer/doxygen/doxfiles/examples/medcouplingexamplesmeshes.doxy @@ -171,7 +171,10 @@ finds one reversed face. After that we fix the incorrectly oriented cell using \ref MEDCoupling::MEDCouplingUMesh::orientCorrectly2DCells "orientCorrectly2DCells()" and re-check the orientation of cells. \snippet MEDCouplingExamplesTest.cxx CppSnippet_MEDCouplingUMesh_are2DCellsNotCorrectlyOriented_2 - +Alternatively you can orient all 2D cells equally using the first cell as a reference: +\snippet MEDCouplingExamplesTest.cxx CppSnippet_MEDCouplingUMesh_are2DCellsNotCorrectlyOriented_3 +Also it is possible to orient some selected 2D cells by using another group of cells as the reference: +\snippet MEDCouplingExamplesTest.cxx CppSnippet_MEDCouplingUMesh_are2DCellsNotCorrectlyOriented_4 \subsubsection cpp_mcumesh_renumberNodesInConn Renumbering nodes in the connectivity array diff --git a/doc/developer/doxygen/fakesources/MEDCouplingUMesh.C b/doc/developer/doxygen/fakesources/MEDCouplingUMesh.C index b9e0ac718..6a4fa25d9 100644 --- a/doc/developer/doxygen/fakesources/MEDCouplingUMesh.C +++ b/doc/developer/doxygen/fakesources/MEDCouplingUMesh.C @@ -132,6 +132,7 @@ MEDCouplingUMesh::isEqualWithoutConsideringStr(const MEDCouplingMesh *other, dou //MEDCouplingUMesh::mergeNodes(double precision, bool& areNodesMerged, int& newNbOfNodes); //MEDCouplingUMesh::mergeNodes2(double precision, bool& areNodesMerged, int& newNbOfNodes); MEDCouplingUMesh::orientCorrectly2DCells(const double *vec, bool polyOnly); +MEDCouplingUMesh::orientCorrectly2DCells(const MEDCouplingUMesh* refFaces = nullptr); MEDCouplingUMesh::orientCorrectlyPolyhedrons(); //MEDCouplingUMesh::renumberNodes(const int *newNodeNumbers, int newNbOfNodes); //MEDCouplingUMesh::renumberNodes2(const int *newNodeNumbers, int newNbOfNodes); diff --git a/doc/user/input/data_analysis.rst b/doc/user/input/data_analysis.rst index ed6753fae..1ac5b2b66 100644 --- a/doc/user/input/data_analysis.rst +++ b/doc/user/input/data_analysis.rst @@ -491,6 +491,18 @@ The last argument if True, only polygons are checked, else, all cells are checke A mesh before applying orientCorrectly2DCells (to the left) and after (to the right) +Alternatively you can orient all 2D cells equally using the first cell as a reference: + +.. literalinclude:: ../../../src/MEDCoupling_Swig/UsersGuideExamplesTest.py + :start-after: UG_CommonHandlingMesh_11_1 + :end-before: UG_CommonHandlingMesh_11_1 + +Also it is possible to orient some selected 2D cells by using another group of cells as the reference: + +.. literalinclude:: ../../../src/MEDCoupling_Swig/UsersGuideExamplesTest.py + :start-after: UG_CommonHandlingMesh_11_2 + :end-before: UG_CommonHandlingMesh_11_2 + If your mesh includes incorrectly oriented polyhedra, the following method can help to fix your mesh: .. literalinclude:: ../../../src/MEDCoupling_Swig/UsersGuideExamplesTest.py diff --git a/src/MEDCoupling/MEDCouplingUMesh.cxx b/src/MEDCoupling/MEDCouplingUMesh.cxx index 0762ba266..a9181f3e7 100755 --- a/src/MEDCoupling/MEDCouplingUMesh.cxx +++ b/src/MEDCoupling/MEDCouplingUMesh.cxx @@ -8527,3 +8527,411 @@ const mcIdType *MEDCouplingUMeshCell::getAllConn(mcIdType& lgth) const else return 0; } + +/// @cond INTERNAL + +namespace MEDCouplingImpl +{ + const mcIdType theUndefID = std::numeric_limits< mcIdType >::max(); //!< undefined cell id + + //================================================================================ + /*! + * \brief Encode a cell id and a mesh index into a code + * \param [in] id - cell id + * \param [in] iMesh - mesh index [0,1] + * \return mcIdType - code + */ + //================================================================================ + + mcIdType encodeID( mcIdType id, int iMesh ) + { + return ( id + 1 ) * ( iMesh ? -1 : 1 ); + } + //================================================================================ + /*! + * \brief Return cell id and mesh index by a given id + * \param [in] id - code of a cell in a mesh + * \param [out] iMesh - returned mesh index + * \return mcIdType - cell id + */ + //================================================================================ + + mcIdType decodeID( mcIdType id, int& iMesh ) + { + iMesh = ( id < 0 ); + return std::abs( id ) - 1; + } + + //================================================================================ + /*! + * \brief return another face sharing two given nodes of a face edge + * \param [in] n0 - 1st node of the edge + * \param [in] n1 - 2nd node of the edge + * \param [in] inputFaceID - face including \a n0 andf \a n2 + * \param [in] mesh - object and reference meshes + * \param [in] revNodal - reverse nodal connectivity of the two meshes + * \param [in] revNodalIndx - index of reverse nodal connectivity of the two meshes + * \param [out] facesByEdge - return another face including \a n0 andf \a n2 + * \param [out] equalFaces - return faces equal to facesByEdge + */ + //================================================================================ + + void getFacesOfEdge( mcIdType n0, + mcIdType n1, + mcIdType inputFaceID, + MEDCouplingUMesh* mesh[], + MCAuto revNodal[], + MCAuto revNodalIndx[], + std::vector< mcIdType >& facesByEdge, + std::vector< mcIdType >& equalFaces) + { + // find faces sharing the both nodes of edge + + facesByEdge.clear(); + size_t prevNbF; // nb faces found in 0-th mesh + for ( int iM = 0; iM < 2; ++iM ) + { + const mcIdType * revInd = revNodalIndx[ iM ]->begin(); + const mcIdType * rev = revNodal [ iM ]->begin(); + + mcIdType nbRevFaces0 = revInd[ n0 + 1 ] - revInd[ n0 ]; + mcIdType nbRevFaces1 = revInd[ n1 + 1 ] - revInd[ n1 ]; + + prevNbF = facesByEdge.size(); + facesByEdge.resize( prevNbF + std::max( nbRevFaces0, nbRevFaces1 )); + + auto it = std::set_intersection( rev + revInd[ n0 ], + rev + revInd[ n0 ] + nbRevFaces0, + rev + revInd[ n1 ], + rev + revInd[ n1 ] + nbRevFaces1, + facesByEdge.begin() + prevNbF ); + facesByEdge.resize( it - facesByEdge.begin() ); + } + + // facesByEdge now contains at least the 'inputFaceID' + // check if there are other faces + + size_t nbF = facesByEdge.size(); + if ( nbF > 1 ) + { + if ( prevNbF > 0 && prevNbF < nbF ) // faces found in both meshes + { + // remove from facesByEdge equal faces in different meshes + const mcIdType *conn [2] = { mesh[0]->getNodalConnectivity()->getConstPointer(), + mesh[1]->getNodalConnectivity()->getConstPointer() }; + const mcIdType *connI[2] = { mesh[0]->getNodalConnectivityIndex()->getConstPointer(), + mesh[1]->getNodalConnectivityIndex()->getConstPointer() }; + for ( size_t i0 = 0; i0 < prevNbF; ++i0 ) + { + if ( facesByEdge[ i0 ] == theUndefID ) + continue; + mcIdType objFaceID = MEDCouplingImpl::encodeID( facesByEdge[ i0 ], 0 ); + bool isInputFace = ( objFaceID == inputFaceID ); + + for ( size_t i1 = prevNbF; i1 < facesByEdge.size(); ++i1 ) + { + if ( facesByEdge[ i1 ] == theUndefID ) + continue; + + mcIdType f0 = facesByEdge[ i0 ]; + mcIdType f1 = facesByEdge[ i1 ]; + size_t nbNodes0 = connI[0][ f0 + 1 ] - connI[0][ f0 ] - 1; + size_t nbNodes1 = connI[1][ f1 + 1 ] - connI[1][ f1 ] - 1; + if ( nbNodes0 != nbNodes1 ) + continue; + + const mcIdType * fConn0 = conn[0] + connI[0][ f0 ] + 1; + const mcIdType * fConn1 = conn[1] + connI[1][ f1 ] + 1; + if ( std::equal( fConn0, fConn0 + nbNodes0, fConn1 )) + { + // equal faces; remove an object one + mcIdType refFaceID = MEDCouplingImpl::encodeID( facesByEdge[ i1 ], 1 ); + if ( refFaceID == inputFaceID ) + isInputFace = true; + + if ( std::find( equalFaces.begin(), + equalFaces.end(), objFaceID ) == equalFaces.end() ) + equalFaces.push_back( objFaceID ); + + facesByEdge[ i0 ] = theUndefID; + if ( isInputFace ) + facesByEdge[ i1 ] = theUndefID; + break; + } + } + if ( isInputFace ) + facesByEdge[ i0 ] = theUndefID; + } + } + } + + nbF = facesByEdge.size(); + for ( size_t i = 0; i < facesByEdge.size(); ++i ) + { + if ( facesByEdge[ i ] != theUndefID ) + { + facesByEdge[ i ] = MEDCouplingImpl::encodeID( facesByEdge[ i ], i >= prevNbF ); + if ( facesByEdge[ i ] == inputFaceID ) + facesByEdge[ i ] = theUndefID; + } + nbF -= ( facesByEdge[ i ] == theUndefID ); + } + + if ( nbF > 1 ) + return; // non-manifold + + if ( nbF < 1 ) + { + facesByEdge.clear(); + } + else // nbF == 1, set a found face first + { + if ( facesByEdge[ 0 ] == theUndefID ) + { + for ( size_t i = 1; i < facesByEdge.size(); ++i ) + if ( facesByEdge[ i ] != theUndefID ) + { + facesByEdge[ 0 ] = facesByEdge[ i ]; + break; + } + } + facesByEdge.resize( 1 ); + } + return; + } + + //================================================================================ + /*! + * \brief Remove a face from nodal reversed connectivity + * \param [in] node - a node of the face + * \param [in] face - the face + * \param [in.out] revNodal - reverse nodal connectivity + * \param [in,out] revNodalIndx - reverse nodal connectivity index + */ + //================================================================================ + + void removeFromRevNodal( mcIdType node, + mcIdType face, + MCAuto& revNodal, + MCAuto& revNodalIndx) + { + mcIdType* fBeg = revNodal->getPointer() + revNodalIndx->getIJ( node, 0 ); + mcIdType* fEnd = revNodal->getPointer() + revNodalIndx->getIJ( node + 1, 0); + auto it = std::find( fBeg, fEnd, face ); + if ( it != fEnd ) + { + for ( auto it2 = it + 1; it2 < fEnd; ++it2 ) // keep faces sorted + *( it2 - 1 ) = *it2; + + *( fEnd - 1 ) = theUndefID; + } + } + + //================================================================================ + /*! + * \brief Check order of two nodes in a given face + * \param [inout] n0 - node 1 + * \param [inout] n1 - node 2 + * \param [inout] iFEnc - face + * \param [inout] mesh - mesh + * \return bool - true if the nodes are in [ .., n1, n0, ..] order in face + */ + //================================================================================ + + bool isReverseOrder( mcIdType n0, + mcIdType n1, + mcIdType iFEnc, + MEDCouplingUMesh* mesh[] ) + { + int iMesh; + mcIdType iF = decodeID( iFEnc, iMesh ); + + const mcIdType *conn = mesh[ iMesh ]->getNodalConnectivity()->getConstPointer(); + const mcIdType *connI = mesh[ iMesh ]->getNodalConnectivityIndex()->getConstPointer(); + + auto it0 = std::find( conn + connI[ iF ] + 1, + conn + connI[ iF + 1 ], + n0 ); + auto it1 = std::find( conn + connI[ iF ] + 1, + conn + connI[ iF + 1 ], + n1 ); + long i0 = it0 - conn; + long i1 = it1 - conn; + + bool isRev = ( std::abs( i1 - i0 ) == 1 ) ? i1 < i0 : i0 < i1; + return isRev; + } + + //================================================================================ + /*! + * \brief Change orientation of a face in one of given meshes + * \param [in] iFEnc - face ID also encoding a mesh index + * \param [in,out] mesh - object and reference meshes + */ + //================================================================================ + + void reverseFace( mcIdType iFEnc, MEDCouplingUMesh* mesh[] ) + { + int iMesh; + mcIdType face = decodeID( iFEnc, iMesh ); + + mcIdType *conn = mesh[ iMesh ]->getNodalConnectivity()->getPointer(); + mcIdType *connI = mesh[ iMesh ]->getNodalConnectivityIndex()->getPointer(); + + const INTERP_KERNEL::CellModel& cm = + INTERP_KERNEL::CellModel::GetCellModel( mesh[iMesh]->getTypeOfCell( face )); + + cm.changeOrientationOf2D( conn + connI[ face ] + 1, + (unsigned int)( connI[ face + 1 ] - connI[ face ] - 1 )); + return; + } +} + +/// @endcond + +//================================================================================ +/*! + * \brief Orient cells of \a this 2D mesh equally to \a refFaces + * \param [in] refFaces - 2D mesh containing correctly oriented faces. It is optional. + * If there are no cells in \a refFaces or it is nullptr, then any face + * in \a this mesh is used as a reference + * \throw If \a this mesh is not well defined. + * \throw If \a this mesh or \refFaces are not 2D. + * \throw If \a this mesh and \refFaces do not share nodes. + * \throw If \a refFaces are not equally oriented. + * \throw If \a this mesh plus \a refFaces together form a non-manifold mesh. + * + * \if ENABLE_EXAMPLES + * \ref cpp_mcumesh_are2DCellsNotCorrectlyOriented "Here is a C++ example".
+ * \ref py_mcumesh_are2DCellsNotCorrectlyOriented "Here is a Python example". + * \endif + */ +//================================================================================ + +void MEDCouplingUMesh::orientCorrectly2DCells(const MEDCouplingUMesh* refFaces) +{ + checkConsistencyLight(); + if ( getMeshDimension() != 2 ) + throw INTERP_KERNEL::Exception("The mesh dimension must be 2"); + if ( refFaces ) + { + refFaces->checkConsistencyLight(); + if ( refFaces->getMeshDimension() != 2 ) + throw INTERP_KERNEL::Exception("The reference mesh dimension must be 2"); + if ( getCoords() != refFaces->getCoords() ) + throw INTERP_KERNEL::Exception("Object and reference meshes must share nodes "); + if ( refFaces->getNumberOfCells() == 0 ) + refFaces = nullptr; + } + if ( getNumberOfCells() == 0 ) + return; + + enum { _OBJ, _REF }; + MEDCouplingUMesh* mesh[2] = { this, const_cast< MEDCouplingUMesh* >( refFaces ) }; + MCAuto meshPtr; + if ( !mesh[_REF] ) + { + meshPtr = mesh[_REF] = MEDCouplingUMesh::New(); + mesh[_REF]->setCoords( mesh[_OBJ]->getCoords() ); + mesh[_REF]->allocateCells(0); + mesh[_REF]->finishInsertingCells(); + } + mcIdType nbFacesToCheck[2] = { mesh[_OBJ]->getNumberOfCells(), + mesh[_REF]->getNumberOfCells() }; + std::vector< bool > isFaceQueued[ 2 ]; // enqueued faces of 2 meshes + isFaceQueued[_OBJ].resize( nbFacesToCheck[_OBJ] ); + isFaceQueued[_REF].resize( nbFacesToCheck[_REF] ); + + MCAuto revNodal [2] = { DataArrayIdType::New(), DataArrayIdType::New() }; + MCAuto revNodalIndx[2] = { DataArrayIdType::New(), DataArrayIdType::New() }; + mesh[_OBJ]->getReverseNodalConnectivity( revNodal[_OBJ], revNodalIndx[_OBJ] ); + mesh[_REF]->getReverseNodalConnectivity( revNodal[_REF], revNodalIndx[_REF] ); + + std::vector< mcIdType > faceNodes(4); + std::vector< mcIdType > facesByEdge(4), equalFaces; + std::vector< mcIdType > faceQueue; // starting faces with IDs counted from 1; negative ID mean a face in ref mesh + + while ( nbFacesToCheck[_OBJ] + nbFacesToCheck[_REF] > 0 ) // until all faces checked + { + if ( faceQueue.empty() ) // all neighbors checked, find more faces to check + { + for ( int iMesh = 1; iMesh >= 0; --iMesh ) // on [ _REF, _OBJ ] + if ( nbFacesToCheck[iMesh] > 0 ) + for ( mcIdType f = 0, nbF = mesh[iMesh]->getNumberOfCells(); f < nbF; ++f ) + if ( !isFaceQueued[iMesh][f] ) + { + faceQueue.push_back( MEDCouplingImpl::encodeID( f, iMesh )); + isFaceQueued[ iMesh ][ f ] = true; + iMesh = 0; + break; + } + if ( faceQueue.empty() ) + break; + } + + mcIdType fID = faceQueue.back(); + faceQueue.pop_back(); + + int iMesh, iMesh2; + mcIdType refFace = MEDCouplingImpl::decodeID( fID, iMesh ); + + nbFacesToCheck[iMesh]--; + + equalFaces.clear(); + faceNodes.clear(); + mesh[iMesh]->getNodeIdsOfCell( refFace, faceNodes ); + const INTERP_KERNEL::CellModel& cm = INTERP_KERNEL::CellModel::GetCellModel( mesh[iMesh]->getTypeOfCell( refFace )); + const int nbEdges = cm.getNumberOfSons(); + + // loop on edges of the refFace + mcIdType n0 = faceNodes[ nbEdges - 1 ]; // 1st node of edge + for ( int edge = 0; edge < nbEdges; ++edge ) + { + mcIdType n1 = faceNodes[ edge ]; // 2nd node of edge + + // get faces sharing the edge + MEDCouplingImpl::getFacesOfEdge( n0, n1, fID, mesh, revNodal, revNodalIndx, + facesByEdge, equalFaces ); + + if ( facesByEdge.size() > 1 ) + THROW_IK_EXCEPTION("Non-manifold mesh at edge " << n0+1 << " - " << n1+1); + + if ( facesByEdge.size() == 1 ) + { + // compare orientation of two faces + // + if ( !MEDCouplingImpl::isReverseOrder( n0, n1, facesByEdge[0], mesh )) + { + if ( facesByEdge[0] < 0 ) // in the ref mesh + throw INTERP_KERNEL::Exception("Different orientation of reference faces"); + + MEDCouplingImpl::reverseFace( facesByEdge[0], mesh ); + } + mcIdType face2 = MEDCouplingImpl::decodeID( facesByEdge[0], iMesh2 ); + if ( !isFaceQueued[iMesh2][face2] ) + { + isFaceQueued[iMesh2][face2] = true; + faceQueue.push_back( facesByEdge[0] ); + } + } + n0 = n1; + } + + // remove face and equalFaces from revNodal in order not to treat them again + equalFaces.push_back( fID ); + for ( mcIdType face : equalFaces ) + { + mcIdType f = MEDCouplingImpl::decodeID( face, iMesh2 ); + const mcIdType *conn = mesh[iMesh2]->getNodalConnectivity()->getConstPointer(); + const mcIdType *connI = mesh[iMesh2]->getNodalConnectivityIndex()->getConstPointer(); + mcIdType nbNodes = connI[ f + 1 ] - connI[ f ] - 1; + for ( const mcIdType* n = conn + connI[ f ] + 1, *nEnd = n + nbNodes; n < nEnd; ++n ) + + MEDCouplingImpl::removeFromRevNodal( *n, f, // not to treat f again + revNodal[ iMesh2 ], revNodalIndx[ iMesh2 ] ); + } + + } // while() until all faces checked + + return; +} diff --git a/src/MEDCoupling/MEDCouplingUMesh.hxx b/src/MEDCoupling/MEDCouplingUMesh.hxx index 2d5f408d7..38989739e 100644 --- a/src/MEDCoupling/MEDCouplingUMesh.hxx +++ b/src/MEDCoupling/MEDCouplingUMesh.hxx @@ -197,6 +197,7 @@ namespace MEDCoupling MEDCOUPLING_EXPORT bool removeDegenerated1DCells(); MEDCOUPLING_EXPORT void are2DCellsNotCorrectlyOriented(const double *vec, bool polyOnly, std::vector& cells) const; MEDCOUPLING_EXPORT void orientCorrectly2DCells(const double *vec, bool polyOnly); + MEDCOUPLING_EXPORT void orientCorrectly2DCells(const MEDCouplingUMesh* refFaces = nullptr); MEDCOUPLING_EXPORT void changeOrientationOfCells(); MEDCOUPLING_EXPORT void arePolyhedronsNotCorrectlyOriented(std::vector& cells) const; MEDCOUPLING_EXPORT void orientCorrectlyPolyhedrons(); diff --git a/src/MEDCoupling/Test/MEDCouplingExamplesTest.cxx b/src/MEDCoupling/Test/MEDCouplingExamplesTest.cxx index 37f145cff..7976969be 100644 --- a/src/MEDCoupling/Test/MEDCouplingExamplesTest.cxx +++ b/src/MEDCoupling/Test/MEDCouplingExamplesTest.cxx @@ -919,6 +919,17 @@ void CppExample_MEDCouplingUMesh_are2DCellsNotCorrectlyOriented() mesh->are2DCellsNotCorrectlyOriented( vec, false, badCellIds ); CPPUNIT_ASSERT( badCellIds.size() == 0 ); // the orientation is OK //! [CppSnippet_MEDCouplingUMesh_are2DCellsNotCorrectlyOriented_2] + //! [CppSnippet_MEDCouplingUMesh_are2DCellsNotCorrectlyOriented_3] + mesh->orientCorrectly2DCells(); + //! [CppSnippet_MEDCouplingUMesh_are2DCellsNotCorrectlyOriented_3] + //! [CppSnippet_MEDCouplingUMesh_are2DCellsNotCorrectlyOriented_4] + const mcIdType refCells[] = { 0,2 }; + const mcIdType objCells[] = { 1,3 }; + MCAuto refGroup = mesh->buildPartOfMySelf( refCells, refCells + 2 ); + MCAuto objGroup = mesh->buildPartOfMySelf( objCells, objCells + 2 ); + objGroup->orientCorrectly2DCells( refGroup ); + mesh->setPartOfMySelf( objCells, objCells + 2, *objGroup ); + //! [CppSnippet_MEDCouplingUMesh_are2DCellsNotCorrectlyOriented_4] } void CppExample_MEDCouplingUMesh_getCellsContainingPoints() diff --git a/src/MEDCoupling_Swig/MEDCouplingBasicsTest1.py b/src/MEDCoupling_Swig/MEDCouplingBasicsTest1.py index 7182bc9bc..b347dbcb9 100644 --- a/src/MEDCoupling_Swig/MEDCouplingBasicsTest1.py +++ b/src/MEDCoupling_Swig/MEDCouplingBasicsTest1.py @@ -2311,6 +2311,69 @@ class MEDCouplingBasicsTest1(unittest.TestCase): self.assertEqual(mesh.getNodalConnectivityIndex().getValues(), cIRef) pass + def testCellOrientation6(self): + # CCTP 2.3.1 (bos #26452) + mesh = MEDCouplingUMesh('Orientation', 3) + coo = [0,0,0, 0,10,0, 10,0,0, 10,10,0, -10,0,0, -10,10,-0, 0,0,10, 10,0,10, -10,0,10, 0,-10,0, 10,-10,0, -10,-10,0, 0,5,0, 5,0,0, 10,5,0, 5,10,0, -5,0,-0, -10,5,-0, -5,10,-0, 0,0,5, 10,0,5, 5,0,10, -10,0,5, -5,0,10, 0,-5,0, 10,-5,0, 5,-10,0, -10,-5,0, -5,-10,0, 5,5,0, -5,5,-0, 5,0,5, -5,0,5, 5,-5,0, -5,-5,0] + meshCoords = DataArrayDouble.New(coo, 35, 3) + mesh.setCoords(meshCoords) + mesh.setMeshDimension( 2 ) + conn = [ 1,12,29,15, 15,29,14, 3, 12, 0,13,29, 29,13, 2,14, 1,12,30,18, 18,30,17, 5, 12, 0,16,30, 30,16, 4,17, 6,19,31,21, 21,31,20, 7, 19, 0,13,31, 31,13, 2,20, 6,19,32,23, 23,32,22, 8, 19, 0,16,32, 32,16, 4,22, 9,24,33,26, 26,33,25,10, 24, 0,13,33, 33,13, 2,25, 9,24,34,28, 28,34,27,11, 24, 0,16,34, 34,16, 4,27] + mesh.allocateCells(24) + for i in range(24): + mesh.insertNextCell(NORM_QUAD4,4,conn[4*i:4*(i+1)]); + mesh.finishInsertingCells() + Group_1 = list( range( 0, 4 )) + Group_2 = list( range( 4, 8 )) + Group_3 = list( range( 8, 12 )) + Group_5 = list( range( 16, 20 )) + Group_7 = Group_3 + Group_8 = list( range( 12, 16 )) + + # example 1.1 - check failure on non-manifold + objMesh = mesh.buildPartOfMySelf( Group_1 + Group_3 + Group_5 ) + refMesh = mesh.buildPartOfMySelf( Group_1 ) + self.assertRaisesRegex(Exception, "Non-manifold", + objMesh.orientCorrectly2DCells, refMesh ) + # example 1.2 + # - do nothing, as request "de non-appartenance de la reference a la cible" dropped + # + # example 1.3 - fix orientation of Group_1 and Group_8 + objMesh = mesh.buildPartOfMySelf( Group_1 + Group_2 + Group_7 + Group_8 ) + refMesh = mesh.buildPartOfMySelf( Group_7 ) + objMesh.orientCorrectly2DCells( refMesh ) + # check Group_1 + self.assertEqual( objMesh.getNodeIdsOfCell( 0 ), [ 1, 15, 29, 12 ]) + self.assertEqual( objMesh.getNodeIdsOfCell( 1 ), [ 15, 3, 14, 29 ]) + self.assertEqual( objMesh.getNodeIdsOfCell( 2 ), [ 12, 29, 13, 0 ]) + self.assertEqual( objMesh.getNodeIdsOfCell( 3 ), [ 29, 14, 2, 13 ]) + # check Group_8 + self.assertEqual( objMesh.getNodeIdsOfCell( 12 ), [ 6, 23, 32, 19 ]) + self.assertEqual( objMesh.getNodeIdsOfCell( 13 ), [ 23, 8, 22, 32 ]) + self.assertEqual( objMesh.getNodeIdsOfCell( 14 ), [ 19, 32, 16, 0 ]) + self.assertEqual( objMesh.getNodeIdsOfCell( 15 ), [ 32, 22, 4, 16 ]) + # + # Case with no reference given. Group_2 and Group_7 must reverse + objMesh = mesh.buildPartOfMySelf( Group_1 + Group_2 + Group_7 + Group_8 ) + objMesh.orientCorrectly2DCells( None ) + # check Group_2 + self.assertEqual( objMesh.getNodeIdsOfCell( 4 ), [ 1, 18, 30, 12 ]) + self.assertEqual( objMesh.getNodeIdsOfCell( 5 ), [ 18, 5, 17, 30 ]) + self.assertEqual( objMesh.getNodeIdsOfCell( 6 ), [ 12, 30, 16, 0 ]) + self.assertEqual( objMesh.getNodeIdsOfCell( 7 ), [ 30, 17, 4, 16 ]) + # check Group_7 + self.assertEqual( objMesh.getNodeIdsOfCell( 8 ), [ 6, 21, 31, 19 ]) + self.assertEqual( objMesh.getNodeIdsOfCell( 9 ), [ 21, 7, 20, 31 ]) + self.assertEqual( objMesh.getNodeIdsOfCell( 10 ), [ 19, 31, 13, 0 ]) + self.assertEqual( objMesh.getNodeIdsOfCell( 11 ), [ 31, 20, 2, 13 ]) + # + # Case with differently oriented reference faces. Expect an exception + objMesh = mesh.buildPartOfMySelf( Group_1 + Group_2 + Group_7 + Group_8 ) + refMesh = mesh.buildPartOfMySelf( Group_1 + Group_2 ) + self.assertRaisesRegex(Exception, "Different orientation", + objMesh.orientCorrectly2DCells, refMesh ) + pass + def testPolyhedronBarycenter(self): connN=[0,3,2,1, -1, 4,5,6,7, -1, 0,4,7,3, -1, 3,7,6,2, -1, 2,6,5,1, -1, 1,5,4,0]; coords=[0.,0.,0., 1.,0.,0., 1.,1.,0., 0.,1.,0., 0.,0.,1., 1.,0.,1., 1.,1.,1., 0.,1.,1., 0.5, 0.5, 0.5]; diff --git a/src/MEDCoupling_Swig/MEDCouplingCommon.i b/src/MEDCoupling_Swig/MEDCouplingCommon.i index c2f9afff3..b8e7d08e6 100644 --- a/src/MEDCoupling_Swig/MEDCouplingCommon.i +++ b/src/MEDCoupling_Swig/MEDCouplingCommon.i @@ -2118,6 +2118,7 @@ namespace MEDCoupling DataArrayDouble *getBoundingBoxForBBTree2DQuadratic(double arcDetEps=1e-12) const; DataArrayDouble *getBoundingBoxForBBTree1DQuadratic(double arcDetEps=1e-12) const; void changeOrientationOfCells(); + void orientCorrectly2DCells(const MEDCouplingUMesh *refFaces); DataArrayDouble *computeCellCenterOfMassWithPrecision(double eps); int split2DCells(const DataArrayIdType *desc, const DataArrayIdType *descI, const DataArrayIdType *subNodesInSeg, const DataArrayIdType *subNodesInSegI, const DataArrayIdType *midOpt=0, const DataArrayIdType *midOptI=0); static MEDCouplingUMesh *Build0DMeshFromCoords(DataArrayDouble *da); diff --git a/src/MEDCoupling_Swig/MEDCouplingExamplesTest.py b/src/MEDCoupling_Swig/MEDCouplingExamplesTest.py index bda6078a2..625472e53 100644 --- a/src/MEDCoupling_Swig/MEDCouplingExamplesTest.py +++ b/src/MEDCoupling_Swig/MEDCouplingExamplesTest.py @@ -706,6 +706,17 @@ class MEDCouplingBasicsTest(unittest.TestCase): badCellIds=mesh.are2DCellsNotCorrectlyOriented( vec, False ) assert len( badCellIds ) == 0 # the orientation is OK #! [PySnippet_MEDCouplingUMesh_are2DCellsNotCorrectlyOriented_2] + #! [PySnippet_MEDCouplingUMesh_are2DCellsNotCorrectlyOriented_3] + mesh.orientCorrectly2DCells( None ) + #! [PySnippet_MEDCouplingUMesh_are2DCellsNotCorrectlyOriented_3] + #! [PySnippet_MEDCouplingUMesh_are2DCellsNotCorrectlyOriented_4] + refCells = [ 0,2 ] + objCells = [ 1,3 ] + refGroup = mesh.buildPartOfMySelf( refCells ) + objGroup = mesh.buildPartOfMySelf( objCells ) + objGroup.orientCorrectly2DCells( refGroup ) + mesh.setPartOfMySelf( objCells, objGroup ) + #! [PySnippet_MEDCouplingUMesh_are2DCellsNotCorrectlyOriented_4] return def testExample_MEDCouplingUMesh_getCellsContainingPoints(self): diff --git a/src/MEDCoupling_Swig/UsersGuideExamplesTest.py b/src/MEDCoupling_Swig/UsersGuideExamplesTest.py index 6a4673459..2b7557f34 100755 --- a/src/MEDCoupling_Swig/UsersGuideExamplesTest.py +++ b/src/MEDCoupling_Swig/UsersGuideExamplesTest.py @@ -475,6 +475,19 @@ vec=[0,0,-1] skin.orientCorrectly2DCells(vec,False) #! [UG_CommonHandlingMesh_11] +#! [UG_CommonHandlingMesh_11_1] +skin.orientCorrectly2DCells( None ) +#! [UG_CommonHandlingMesh_11_1] + +#! [UG_CommonHandlingMesh_11_2] +refCells = [ 0,2,4 ] +objCells = [ 1,3,5,6,7,8, 20 ] +refGroup = skin.buildPartOfMySelf( refCells ) +objGroup = skin.buildPartOfMySelf( objCells ) +objGroup.orientCorrectly2DCells( refGroup ) +skin.setPartOfMySelf( objCells, objGroup ) +#! [UG_CommonHandlingMesh_11_2] + #! [UG_CommonHandlingMesh_12] m3.orientCorrectlyPolyhedrons() #! [UG_CommonHandlingMesh_12] -- 2.39.2