X-Git-Url: http://git.salome-platform.org/gitweb/?p=modules%2Fsmesh.git;a=blobdiff_plain;f=src%2FSMESH%2FSMESH_MeshEditor.cxx;h=be874359b73c68b2ab7891c3acb3c09104a2984c;hp=f97e645dc67e8ef3377c411a7affeedc7fefb8ae;hb=499f29d24922cec66e41b41a0039a954993bc6df;hpb=3b4c34ed3671faa92d8554fe2e2b4cee64d9774e diff --git a/src/SMESH/SMESH_MeshEditor.cxx b/src/SMESH/SMESH_MeshEditor.cxx index f97e645dc..be874359b 100644 --- a/src/SMESH/SMESH_MeshEditor.cxx +++ b/src/SMESH/SMESH_MeshEditor.cxx @@ -1,4 +1,4 @@ -// Copyright (C) 2007-2021 CEA/DEN, EDF R&D, OPEN CASCADE +// Copyright (C) 2007-2022 CEA/DEN, EDF R&D, OPEN CASCADE // // Copyright (C) 2003-2007 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN, // CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS @@ -412,8 +412,8 @@ SMDS_MeshElement* SMESH_MeshEditor::AddElement(const vector & nodeIDs, // Modify a compute state of sub-meshes which become empty //======================================================================= -smIdType SMESH_MeshEditor::Remove (const list< smIdType >& theIDs, - const bool isNodes ) +smIdType SMESH_MeshEditor::Remove (const std::list< smIdType >& theIDs, + const bool isNodes ) { ClearLastCreated(); @@ -474,10 +474,87 @@ smIdType SMESH_MeshEditor::Remove (const list< smIdType >& theIDs, return removed; } +//================================================================================ +/*! + * \brief Remove a node and fill a hole appeared, by changing surrounding faces + */ +//================================================================================ + +void SMESH_MeshEditor::RemoveNodeWithReconnection( const SMDS_MeshNode* node ) +{ + if ( ! node ) + return; + + if ( node->NbInverseElements( SMDSAbs_Volume ) > 0 ) + throw SALOME_Exception( "RemoveNodeWithReconnection() applies to 2D mesh only" ); + + // check that only triangles surround the node + for ( SMDS_ElemIteratorPtr fIt = node->GetInverseElementIterator( SMDSAbs_Face ); fIt->more(); ) + { + const SMDS_MeshElement* face = fIt->next(); + if ( face->GetGeomType() != SMDSGeom_TRIANGLE ) + throw SALOME_Exception( "RemoveNodeWithReconnection() applies to triangle mesh only" ); + if ( face->IsQuadratic() ) + throw SALOME_Exception( "RemoveNodeWithReconnection() applies to linear mesh only" ); + } + + std::vector< const SMDS_MeshNode*> neighbours(2); + SMESH_MeshAlgos::IsOn2DBoundary( node, & neighbours ); + + bool toRemove = ( neighbours.size() > 2 ); // non-manifold ==> just remove + + // if ( neighbours.size() == 2 ) // on boundary + // { + // // check if theNode and neighbours are on a line + // gp_Pnt pN = SMESH_NodeXYZ( node ); + // gp_Pnt p0 = SMESH_NodeXYZ( neighbours[0] ); + // gp_Pnt p1 = SMESH_NodeXYZ( neighbours[1] ); + // double dist01 = p0.Distance( p1 ); + // double tol = 0.01 * dist01; + // double distN = ( gp_Vec( p0, p1 ) ^ gp_Vec( p0, pN )).Magnitude() / dist01; + // bool onLine = distN < tol; + // toRemove = !onLine; + // } + + if ( neighbours.empty() ) // not on boundary + { + TIDSortedElemSet linkedNodes; + GetLinkedNodes( node, linkedNodes, SMDSAbs_Face ); + for ( const SMDS_MeshElement* e : linkedNodes ) neighbours.push_back( cast2Node( e )); + if ( neighbours.empty() ) + toRemove = true; + } + + if ( toRemove ) + { + this->Remove( std::list< smIdType >( 1, node->GetID() ), /*isNode=*/true ); + return; + } + + // choose a node to replace by + const SMDS_MeshNode* nToReplace = nullptr; + SMESH_NodeXYZ nodeXYZ = node; + double minDist = Precision::Infinite(); + for ( const SMDS_MeshNode* n : neighbours ) + { + double dist = nodeXYZ.SquareDistance( n ); + if ( dist < minDist ) + { + minDist = dist; + nToReplace = n; + } + } + + // remove node + replace by nToReplace + std::list< const SMDS_MeshNode* > nodeGroup = { nToReplace, node }; + TListOfListOfNodes nodesToMerge( 1, nodeGroup ); + this->MergeNodes( nodesToMerge ); +} + //================================================================================ /*! * \brief Create 0D elements on all nodes of the given object. - * \param elements - Elements on whose nodes to create 0D elements; if empty, + * \param elements - Elements on whose nodes to create 0D elements; if empty, * the all mesh is treated * \param all0DElems - returns all 0D elements found or created on nodes of \a elements * \param duplicateElements - to add one more 0D element to a node or not @@ -947,24 +1024,24 @@ bool getQuadrangleNodes(const SMDS_MeshNode * theQuadNodes [], { if( tr1->NbNodes() != tr2->NbNodes() ) return false; + // find the 4-th node to insert into tr1 const SMDS_MeshNode* n4 = 0; SMDS_ElemIteratorPtr it = tr2->nodesIterator(); - int i=0; - while ( !n4 && i<3 ) { + for ( int i = 0; !n4 && i < 3; ++i ) + { const SMDS_MeshNode * n = cast2Node( it->next() ); - i++; bool isDiag = ( n == theNode1 || n == theNode2 ); if ( !isDiag ) n4 = n; } + // Make an array of nodes to be in a quadrangle int iNode = 0, iFirstDiag = -1; it = tr1->nodesIterator(); - i=0; - while ( i<3 ) { + for ( int i = 0; i < 3; ++i ) + { const SMDS_MeshNode * n = cast2Node( it->next() ); - i++; bool isDiag = ( n == theNode1 || n == theNode2 ); if ( isDiag ) { if ( iFirstDiag < 0 ) @@ -1079,6 +1156,210 @@ bool SMESH_MeshEditor::DeleteDiag (const SMDS_MeshNode * theNode1, return true; } +//======================================================================= +//function : SplitEdge +//purpose : Replace each triangle bound by theNode1-theNode2 segment with +// two triangles by connecting a node made on the link with a node opposite to the link. +//======================================================================= + +void SMESH_MeshEditor::SplitEdge (const SMDS_MeshNode * theNode1, + const SMDS_MeshNode * theNode2, + double thePosition) +{ + ClearLastCreated(); + + SMESHDS_Mesh * mesh = GetMeshDS(); + + // Get triangles and segments to divide + + std::vector diagNodes = { theNode1, theNode2 }; + std::vector foundElems; + if ( !mesh->GetElementsByNodes( diagNodes, foundElems ) || foundElems.empty() ) + throw SALOME_Exception( SMESH_Comment("No triangle is bound by the edge ") + << theNode1->GetID() << " - " << theNode2->GetID()); + + SMESH_MesherHelper helper( *GetMesh() ); + + for ( const SMDS_MeshElement * elem : foundElems ) + { + SMDSAbs_ElementType type = elem->GetType(); + switch ( type ) { + case SMDSAbs_Volume: + throw SALOME_Exception( "Can't split an edge of a volume"); + break; + + case SMDSAbs_Face: + if ( elem->GetGeomType() != SMDSGeom_TRIANGLE ) + throw SALOME_Exception( "Can't split an edge of a face of type other than triangle"); + if ( elem->IsQuadratic() ) + { + helper.SetIsQuadratic( true ); + helper.AddTLinks( static_cast< const SMDS_MeshFace*>( elem )); + helper.SetIsBiQuadratic( elem->GetEntityType() == SMDSEntity_BiQuad_Triangle ); + } + break; + + case SMDSAbs_Edge: + if ( elem->IsQuadratic() ) + { + helper.SetIsQuadratic( true ); + helper.AddTLinks( static_cast< const SMDS_MeshEdge*>( elem )); + } + break; + default:; + } + } + + // Make a new node + + const SMDS_MeshNode* nodeOnLink = helper.GetMediumNode( theNode1, theNode2,/*force3d=*/false ); + + gp_Pnt newNodeXYZ = ( SMESH_NodeXYZ( theNode1 ) * ( 1 - thePosition ) + + SMESH_NodeXYZ( theNode2 ) * thePosition ); + + const TopoDS_Shape& S = mesh->IndexToShape( nodeOnLink->GetShapeID() ); + if ( !S.IsNull() && S.ShapeType() == TopAbs_FACE ) // find newNodeXYZ by UV on FACE + { + Handle(ShapeAnalysis_Surface) surface = helper.GetSurface( TopoDS::Face( S )); + double tol = 100 * helper.MaxTolerance( S ); + gp_Pnt2d uv = surface->ValueOfUV( newNodeXYZ, tol ); + if ( surface->Gap() < SMESH_NodeXYZ( theNode1 ).Distance( theNode2 )) + { + newNodeXYZ = surface->Value( uv ); + if ( SMDS_FacePositionPtr nPos = nodeOnLink->GetPosition()) + nPos->SetParameters( uv.X(), uv.Y() ); + } + } + if ( !S.IsNull() && S.ShapeType() == TopAbs_EDGE ) // find newNodeXYZ by param on EDGE + { + mesh->MoveNode( nodeOnLink, newNodeXYZ.X(), newNodeXYZ.Y(), newNodeXYZ.Z() ); + double u = Precision::Infinite(), tol = 100 * helper.MaxTolerance( S ), distXYZ[4]; + helper.ToFixNodeParameters( true ); + if ( helper.CheckNodeU( TopoDS::Edge( S ), nodeOnLink, u, tol, /*force3D=*/false, distXYZ )) + newNodeXYZ.SetCoord( distXYZ[1], distXYZ[2], distXYZ[3] ); + } + mesh->MoveNode( nodeOnLink, newNodeXYZ.X(), newNodeXYZ.Y(), newNodeXYZ.Z() ); + + // Split triangles and segments + + std::vector nodes( 7 ); + for ( const SMDS_MeshElement * elem : foundElems ) + { + nodes.assign( elem->begin_nodes(), elem->end_nodes() ); + nodes.resize( elem->NbCornerNodes() + 1 ); + nodes.back() = nodes[0]; + + smIdType id = elem->GetID(); + int shapeID = elem->GetShapeID(); + + const SMDS_MeshNode* centralNode = nullptr; + if ( elem->GetEntityType() == SMDSEntity_BiQuad_Triangle ) + centralNode = elem->GetNode( 6 ); + + mesh->RemoveFreeElement( elem, /*sm=*/0, /*fromGroups=*/false ); + if ( centralNode ) + mesh->RemoveFreeNode( centralNode, /*sm=*/0, /*fromGroups=*/true ); + + for ( size_t i = 1; i < nodes.size(); ++i ) + { + const SMDS_MeshNode* n1 = nodes[i-1]; + const SMDS_MeshNode* n2 = nodes[i]; + const SMDS_MeshElement* newElem; + if ( nodes.size() == 4 ) // triangle + { + bool isDiag1 = ( n1 == theNode1 || n1 == theNode2 ); + bool isDiag2 = ( n2 == theNode1 || n2 == theNode2 ); + if ( isDiag1 && isDiag2 ) + continue; + + newElem = helper.AddFace( n1, n2, nodeOnLink, id ); + } + else // segment + { + newElem = helper.AddEdge( n1, nodeOnLink, id ); + } + myLastCreatedElems.push_back( newElem ); + AddToSameGroups( newElem, elem, mesh ); + if ( shapeID ) + mesh->SetMeshElementOnShape( newElem, shapeID ); + id = 0; + } + } + return; +} + +//======================================================================= +//function : SplitFace +//purpose : Split a face into triangles each formed by two nodes of the +// face and a new node added at the given coordinates. +//======================================================================= + +void SMESH_MeshEditor::SplitFace (const SMDS_MeshElement * theFace, + double theX, + double theY, + double theZ ) +{ + ClearLastCreated(); + + if ( !theFace ) + throw SALOME_Exception("Null face given"); + if ( theFace->GetType() != SMDSAbs_Face ) + throw SALOME_Exception("Not a face given"); + + SMESHDS_Mesh * mesh = GetMeshDS(); + + SMESH_MesherHelper helper( *GetMesh() ); + if ( theFace->IsQuadratic() ) + { + helper.SetIsQuadratic( true ); + helper.AddTLinks( static_cast< const SMDS_MeshFace*>( theFace )); + } + const TopoDS_Shape& shape = mesh->IndexToShape( theFace->GetShapeID() ); + helper.SetSubShape( shape ); + helper.SetElementsOnShape( true ); + + // Make a new node + + const SMDS_MeshNode* centralNode = nullptr; + if ( theFace->GetEntityType() == SMDSEntity_BiQuad_Triangle ) + centralNode = theFace->GetNode( 6 ); + else if ( theFace->GetEntityType() == SMDSEntity_BiQuad_Quadrangle ) + centralNode = theFace->GetNode( 8 ); + + if ( centralNode ) + { + helper.SetIsBiQuadratic( true ); + mesh->MoveNode( centralNode, theX, theY, theZ ); + } + else + centralNode = helper.AddNode( theX, theY, theZ ); + + + // Split theFace + + std::vector nodes( theFace->NbNodes() + 1 ); + nodes.assign( theFace->begin_nodes(), theFace->end_nodes() ); + nodes.resize( theFace->NbCornerNodes() + 1 ); + nodes.back() = nodes[0]; + + smIdType id = theFace->GetID(); + int shapeID = theFace->GetShapeID(); + + mesh->RemoveFreeElement( theFace, /*sm=*/0, /*fromGroups=*/false ); + + for ( size_t i = 1; i < nodes.size(); ++i ) + { + const SMDS_MeshElement* newElem = helper.AddFace( nodes[i-1], nodes[i], centralNode, id ); + + myLastCreatedElems.push_back( newElem ); + AddToSameGroups( newElem, theFace, mesh ); + if ( shapeID ) + mesh->SetMeshElementOnShape( newElem, shapeID ); + id = 0; + } + return; +} + //======================================================================= //function : Reorient //purpose : Reverse theElement orientation @@ -1161,69 +1442,88 @@ bool SMESH_MeshEditor::Reorient (const SMDS_MeshElement * theElem) /*! * \brief Reorient faces. * \param theFaces - the faces to reorient. If empty the whole mesh is meant - * \param theDirection - desired direction of normal of \a theFace - * \param theFace - one of \a theFaces that should be oriented according to - * \a theDirection and whose orientation defines orientation of other faces + * \param theDirection - desired direction of normal of \a theRefFaces. + * It can be (0,0,0) in order to keep orientation of \a theRefFaces. + * \param theRefFaces - correctly oriented faces whose orientation defines + * orientation of other faces. * \return number of reoriented faces. */ //================================================================================ -int SMESH_MeshEditor::Reorient2D (TIDSortedElemSet & theFaces, - const gp_Dir& theDirection, - const SMDS_MeshElement * theFace) +int SMESH_MeshEditor::Reorient2D( TIDSortedElemSet & theFaces, + const gp_Vec& theDirection, + TIDSortedElemSet & theRefFaces, + bool theAllowNonManifold ) { int nbReori = 0; - if ( !theFace || theFace->GetType() != SMDSAbs_Face ) return nbReori; if ( theFaces.empty() ) { - SMDS_FaceIteratorPtr fIt = GetMeshDS()->facesIterator(/*idInceasingOrder=true*/); + SMDS_FaceIteratorPtr fIt = GetMeshDS()->facesIterator(); while ( fIt->more() ) theFaces.insert( theFaces.end(), fIt->next() ); + + if ( theFaces.empty() ) + return nbReori; } - // orient theFace according to theDirection - gp_XYZ normal; - SMESH_MeshAlgos::FaceNormal( theFace, normal, /*normalized=*/false ); - if ( normal * theDirection.XYZ() < 0 ) - nbReori += Reorient( theFace ); + // orient theRefFaces according to theDirection + if ( theDirection.X() != 0 || theDirection.Y() != 0 || theDirection.Z() != 0 ) + for ( const SMDS_MeshElement* refFace : theRefFaces ) + { + gp_XYZ normal; + SMESH_MeshAlgos::FaceNormal( refFace, normal, /*normalized=*/false ); + if ( normal * theDirection.XYZ() < 0 ) + nbReori += Reorient( refFace ); + } - // Orient other faces + // mark reference faces + GetMeshDS()->SetAllCellsNotMarked(); + for ( const SMDS_MeshElement* refFace : theRefFaces ) + refFace->setIsMarked( true ); - set< const SMDS_MeshElement* > startFaces, visitedFaces; - TIDSortedElemSet avoidSet; - set< SMESH_TLink > checkedLinks; - pair< set< SMESH_TLink >::iterator, bool > linkIt_isNew; + // erase reference faces from theFaces + for ( TIDSortedElemSet::iterator fIt = theFaces.begin(); fIt != theFaces.end(); ) + if ( (*fIt)->isMarked() ) + fIt = theFaces.erase( fIt ); + else + ++fIt; - if ( theFaces.size() > 1 )// leave 1 face to prevent finding not selected faces - theFaces.erase( theFace ); - startFaces.insert( theFace ); + if ( theRefFaces.empty() ) + { + theRefFaces.insert( *theFaces.begin() ); + theFaces.erase( theFaces.begin() ); + } + + // Orient theFaces + + // if ( theFaces.size() > 1 )// leave 1 face to prevent finding not selected faces + // theFaces.erase( theFace ); int nodeInd1, nodeInd2; - const SMDS_MeshElement* otherFace; + const SMDS_MeshElement* refFace, *otherFace; vector< const SMDS_MeshElement* > facesNearLink; vector< std::pair< int, int > > nodeIndsOfFace; + TIDSortedElemSet avoidSet, emptySet; + NCollection_Map< SMESH_TLink, SMESH_TLink > checkedLinks; - set< const SMDS_MeshElement* >::iterator startFace = startFaces.begin(); - while ( !startFaces.empty() ) + while ( !theRefFaces.empty() ) { - startFace = startFaces.begin(); - theFace = *startFace; - startFaces.erase( startFace ); - if ( !visitedFaces.insert( theFace ).second ) - continue; + auto refFaceIt = theRefFaces.begin(); + refFace = *refFaceIt; + theRefFaces.erase( refFaceIt ); avoidSet.clear(); - avoidSet.insert(theFace); + avoidSet.insert( refFace ); - NLink link( theFace->GetNode( 0 ), (SMDS_MeshNode *) 0 ); + NLink link( refFace->GetNode( 0 ), nullptr ); - const int nbNodes = theFace->NbCornerNodes(); - for ( int i = 0; i < nbNodes; ++i ) // loop on links of theFace + const int nbNodes = refFace->NbCornerNodes(); + for ( int i = 0; i < nbNodes; ++i ) // loop on links of refFace { - link.second = theFace->GetNode(( i+1 ) % nbNodes ); - linkIt_isNew = checkedLinks.insert( link ); - if ( !linkIt_isNew.second ) + link.second = refFace->GetNode(( i+1 ) % nbNodes ); + bool isLinkVisited = checkedLinks.Contains( link ); + if ( isLinkVisited ) { // link has already been checked and won't be encountered more // if the group (theFaces) is manifold @@ -1231,28 +1531,41 @@ int SMESH_MeshEditor::Reorient2D (TIDSortedElemSet & theFaces, } else { + checkedLinks.Add( link ); + facesNearLink.clear(); nodeIndsOfFace.clear(); + TIDSortedElemSet::iterator objFaceIt = theFaces.end(); + while (( otherFace = SMESH_MeshAlgos::FindFaceInSet( link.first, link.second, - theFaces, avoidSet, + emptySet, avoidSet, &nodeInd1, &nodeInd2 ))) - if ( otherFace != theFace) + { + if (( otherFace->isMarked() ) || // ref face + (( objFaceIt = theFaces.find( otherFace )) != theFaces.end() )) // object face { facesNearLink.push_back( otherFace ); nodeIndsOfFace.push_back( make_pair( nodeInd1, nodeInd2 )); - avoidSet.insert( otherFace ); } + avoidSet.insert( otherFace ); + } if ( facesNearLink.size() > 1 ) { // NON-MANIFOLD mesh shell ! - // select a face most co-directed with theFace, + if ( !theAllowNonManifold ) + { + throw SALOME_Exception("Non-manifold topology of groups"); + } + // select a face most co-directed with refFace, // other faces won't be visited this time gp_XYZ NF, NOF; - SMESH_MeshAlgos::FaceNormal( theFace, NF, /*normalized=*/false ); + SMESH_MeshAlgos::FaceNormal( refFace, NF, /*normalized=*/false ); double proj, maxProj = -1; - for ( size_t i = 0; i < facesNearLink.size(); ++i ) { + for ( size_t i = 0; i < facesNearLink.size(); ++i ) + { SMESH_MeshAlgos::FaceNormal( facesNearLink[i], NOF, /*normalized=*/false ); - if (( proj = Abs( NF * NOF )) > maxProj ) { + if (( proj = Abs( NF * NOF )) > maxProj ) + { maxProj = proj; otherFace = facesNearLink[i]; nodeInd1 = nodeIndsOfFace[i].first; @@ -1260,9 +1573,9 @@ int SMESH_MeshEditor::Reorient2D (TIDSortedElemSet & theFaces, } } // not to visit rejected faces - for ( size_t i = 0; i < facesNearLink.size(); ++i ) - if ( facesNearLink[i] != otherFace && theFaces.size() > 1 ) - visitedFaces.insert( facesNearLink[i] ); + // for ( size_t i = 0; i < facesNearLink.size(); ++i ) + // if ( facesNearLink[i] != otherFace && theFaces.size() > 1 ) + // visitedFaces.insert( facesNearLink[i] ); } else if ( facesNearLink.size() == 1 ) { @@ -1270,20 +1583,36 @@ int SMESH_MeshEditor::Reorient2D (TIDSortedElemSet & theFaces, nodeInd1 = nodeIndsOfFace.back().first; nodeInd2 = nodeIndsOfFace.back().second; } - if ( otherFace && otherFace != theFace) + if ( otherFace ) { - // link must be reverse in otherFace if orientation to otherFace - // is same as that of theFace - if ( abs(nodeInd2-nodeInd1) == 1 ? nodeInd2 > nodeInd1 : nodeInd1 > nodeInd2 ) + // link must be reverse in otherFace if orientation of otherFace + // is same as that of refFace + if ( abs( nodeInd2 - nodeInd1 ) == 1 ? nodeInd2 > nodeInd1 : nodeInd1 > nodeInd2 ) { + if ( otherFace->isMarked() ) + throw SALOME_Exception("Different orientation of reference faces"); nbReori += Reorient( otherFace ); } - startFaces.insert( otherFace ); + if ( !otherFace->isMarked() ) + { + theRefFaces.insert( otherFace ); + if ( objFaceIt != theFaces.end() ) + theFaces.erase( objFaceIt ); + } } } - std::swap( link.first, link.second ); // reverse the link + link.first = link.second; // reverse the link + + } // loop on links of refFace + + if ( theRefFaces.empty() && !theFaces.empty() ) + { + theRefFaces.insert( *theFaces.begin() ); + theFaces.erase( theFaces.begin() ); } - } + + } // while ( !theRefFaces.empty() ) + return nbReori; }