X-Git-Url: http://git.salome-platform.org/gitweb/?p=modules%2Fsmesh.git;a=blobdiff_plain;f=src%2FSMESH%2FSMESH_MeshEditor.cxx;h=22695bac0567cc812b241097973dc6ae7a5bc028;hp=ceaeb1acb681d86934f3c211e4b07704a6a7622b;hb=HEAD;hpb=00bd895f4432dd67f39c9f4e5066bfdc2f337a24 diff --git a/src/SMESH/SMESH_MeshEditor.cxx b/src/SMESH/SMESH_MeshEditor.cxx index ceaeb1acb..125bbaa2e 100644 --- a/src/SMESH/SMESH_MeshEditor.cxx +++ b/src/SMESH/SMESH_MeshEditor.cxx @@ -1,4 +1,4 @@ -// Copyright (C) 2007-2016 CEA/DEN, EDF R&D, OPEN CASCADE +// Copyright (C) 2007-2024 CEA, EDF, OPEN CASCADE // // Copyright (C) 2003-2007 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN, // CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS @@ -62,6 +62,8 @@ #include #include #include +#include +#include #include #include #include @@ -97,10 +99,12 @@ #include #include -#include #include "SMESH_TryCatch.hxx" // include after OCCT headers! +#include +#include + #define cast2Node(elem) static_cast( elem ) using namespace std; @@ -162,8 +166,7 @@ SMESH_MeshEditor::ElemFeatures::Init( const SMDS_MeshElement* elem, bool basicOn myIsQuad = elem->IsQuadratic(); if ( myType == SMDSAbs_Volume && !basicOnly ) { - vector quant = static_cast( elem )->GetQuantities(); - myPolyhedQuantities.swap( quant ); + myPolyhedQuantities = static_cast( elem )->GetQuantities(); } } } @@ -188,7 +191,7 @@ SMESH_MeshEditor::AddElement(const vector & node, SMDS_MeshElement* e = 0; int nbnode = node.size(); SMESHDS_Mesh* mesh = GetMeshDS(); - const int ID = features.myID; + const smIdType ID = features.myID; switch ( features.myType ) { case SMDSAbs_Face: @@ -298,6 +301,18 @@ SMESH_MeshEditor::AddElement(const vector & node, node[8], node[9], node[10],node[11], node[12],node[13],node[14] ); } + else if (nbnode == 18) { + if ( ID >= 1 ) e = mesh->AddVolumeWithID(node[0], node[1], node[2], node[3], + node[4], node[5], node[6], node[7], + node[8], node[9], node[10],node[11], + node[12],node[13],node[14], + node[15],node[16],node[17],ID ); + else e = mesh->AddVolume (node[0], node[1], node[2], node[3], + node[4], node[5], node[6], node[7], + node[8], node[9], node[10],node[11], + node[12],node[13],node[14], + node[15],node[16],node[17] ); + } else if (nbnode == 20) { if ( ID >= 1 ) e = mesh->AddVolumeWithID(node[0], node[1], node[2], node[3], node[4], node[5], node[6], node[7], @@ -379,12 +394,12 @@ SMESH_MeshEditor::AddElement(const vector & node, */ //======================================================================= -SMDS_MeshElement* SMESH_MeshEditor::AddElement(const vector & nodeIDs, - const ElemFeatures& features) +SMDS_MeshElement* SMESH_MeshEditor::AddElement(const vector & nodeIDs, + const ElemFeatures& features) { vector nodes; nodes.reserve( nodeIDs.size() ); - vector::const_iterator id = nodeIDs.begin(); + vector::const_iterator id = nodeIDs.begin(); while ( id != nodeIDs.end() ) { if ( const SMDS_MeshNode* node = GetMeshDS()->FindNode( *id++ )) nodes.push_back( node ); @@ -400,16 +415,16 @@ SMDS_MeshElement* SMESH_MeshEditor::AddElement(const vector & nodeIDs, // Modify a compute state of sub-meshes which become empty //======================================================================= -int SMESH_MeshEditor::Remove (const list< int >& theIDs, - const bool isNodes ) +smIdType SMESH_MeshEditor::Remove (const std::list< smIdType >& theIDs, + const bool isNodes ) { ClearLastCreated(); SMESHDS_Mesh* aMesh = GetMeshDS(); set< SMESH_subMesh *> smmap; - int removed = 0; - list::const_iterator it = theIDs.begin(); + smIdType removed = 0; + list::const_iterator it = theIDs.begin(); for ( ; it != theIDs.end(); it++ ) { const SMDS_MeshElement * elem; if ( isNodes ) @@ -462,10 +477,87 @@ int SMESH_MeshEditor::Remove (const list< int >& 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 @@ -607,25 +699,6 @@ static void shiftNodesQuadTria(vector< const SMDS_MeshNode* >& aNodes) aNodes[5] = nd2; } -//======================================================================= -//function : nbEdgeConnectivity -//purpose : return number of the edges connected with the theNode. -// if theEdges has connections with the other type of the -// elements, return -1 -//======================================================================= - -static int nbEdgeConnectivity(const SMDS_MeshNode* theNode) -{ - // SMDS_ElemIteratorPtr elemIt = theNode->GetInverseElementIterator(); - // int nb=0; - // while(elemIt->more()) { - // elemIt->next(); - // nb++; - // } - // return nb; - return theNode->NbInverseElements(); -} - //======================================================================= //function : getNodesFromTwoTria //purpose : @@ -954,24 +1027,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 ) @@ -1086,6 +1159,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 @@ -1113,19 +1390,37 @@ bool SMESH_MeshEditor::Reorient (const SMDS_MeshElement * theElem) MESSAGE("Warning: bad volumic element"); return false; } - const int nbFaces = aPolyedre->NbFaces(); + SMDS_VolumeTool vTool( aPolyedre ); + const int nbFaces = vTool.NbFaces(); + vector quantities( nbFaces ); vector poly_nodes; - vector quantities (nbFaces); - // reverse each face of the polyedre - for (int iface = 1; iface <= nbFaces; iface++) { - int inode, nbFaceNodes = aPolyedre->NbFaceNodes(iface); - quantities[iface - 1] = nbFaceNodes; + // check if all facets are oriented equally + bool sameOri = true; + vector& facetOri = quantities; // keep orientation in quantities so far + for (int iface = 0; iface < nbFaces; iface++) + { + facetOri[ iface ] = vTool.IsFaceExternal( iface ); + if ( facetOri[ iface ] != facetOri[ 0 ]) + sameOri = false; + } + + // reverse faces of the polyhedron + int neededOri = sameOri ? 1 - facetOri[0] : 1; + poly_nodes.reserve( vTool.NbNodes() ); + for ( int iface = 0; iface < nbFaces; iface++ ) + { + int nbFaceNodes = vTool.NbFaceNodes( iface ); + const SMDS_MeshNode** nodes = vTool.GetFaceNodes( iface ); + bool toReverse = ( facetOri[ iface ] != neededOri ); - for (inode = nbFaceNodes; inode >= 1; inode--) { - const SMDS_MeshNode* curNode = aPolyedre->GetFaceNode(iface, inode); - poly_nodes.push_back(curNode); - } + quantities[ iface ] = nbFaceNodes; + + if ( toReverse ) + for ( int inode = nbFaceNodes - 1; inode >= 0; inode-- ) + poly_nodes.push_back( nodes[ inode ]); + else + poly_nodes.insert( poly_nodes.end(), nodes, nodes + nbFaceNodes ); } return GetMeshDS()->ChangePolyhedronNodes( theElem, poly_nodes, quantities ); } @@ -1150,69 +1445,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_TLinkHasher > 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 @@ -1220,28 +1534,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; @@ -1249,9 +1576,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 ) { @@ -1259,20 +1586,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; } @@ -1495,9 +1838,12 @@ void SMESH_MeshEditor::QuadTo4Tri (TIDSortedElemSet & theElems) SMESH_MesherHelper helper( *GetMesh() ); helper.SetElementsOnShape( true ); - SMDS_ElemIteratorPtr faceIt; - if ( theElems.empty() ) faceIt = GetMeshDS()->elementsIterator(SMDSAbs_Face); - else faceIt = SMESHUtils::elemSetIterator( theElems ); + // get standalone groups of faces + vector< SMDS_MeshGroup* > allFaceGroups, faceGroups; + for ( SMESHDS_GroupBase* grBase : GetMeshDS()->GetGroups() ) + if ( SMESHDS_Group* group = dynamic_cast( grBase )) + if ( group->GetType() == SMDSAbs_Face && !group->IsEmpty() ) + allFaceGroups.push_back( & group->SMDSGroup() ); bool checkUV; gp_XY uv [9]; uv[8] = gp_XY(0,0); @@ -1508,6 +1854,10 @@ void SMESH_MeshEditor::QuadTo4Tri (TIDSortedElemSet & theElems) Handle(Geom_Surface) surface; TopLoc_Location loc; + SMDS_ElemIteratorPtr faceIt; + if ( theElems.empty() ) faceIt = GetMeshDS()->elementsIterator(SMDSAbs_Face); + else faceIt = SMESHUtils::elemSetIterator( theElems ); + while ( faceIt->more() ) { const SMDS_MeshElement* quad = faceIt->next(); @@ -1584,13 +1934,19 @@ void SMESH_MeshEditor::QuadTo4Tri (TIDSortedElemSet & theElems) myLastCreatedNodes.push_back( nCentral ); } - // create 4 triangles - helper.SetIsQuadratic ( nodes.size() > 4 ); helper.SetIsBiQuadratic( nodes.size() == 9 ); if ( helper.GetIsQuadratic() ) helper.AddTLinks( static_cast< const SMDS_MeshFace*>( quad )); + // select groups to update + faceGroups.clear(); + for ( SMDS_MeshGroup* group : allFaceGroups ) + if ( group->Remove( quad )) + faceGroups.push_back( group ); + + // create 4 triangles + GetMeshDS()->RemoveFreeElement( quad, subMeshDS, /*fromGroups=*/false ); for ( int i = 0; i < 4; ++i ) @@ -1598,8 +1954,9 @@ void SMESH_MeshEditor::QuadTo4Tri (TIDSortedElemSet & theElems) SMDS_MeshElement* tria = helper.AddFace( nodes[ i ], nodes[(i+1)%4], nCentral ); - ReplaceElemInGroups( tria, quad, GetMeshDS() ); myLastCreatedElems.push_back( tria ); + for ( SMDS_MeshGroup* group : faceGroups ) + group->Add( tria ); } } } @@ -1785,6 +2142,17 @@ namespace TSplitMethod( int nbTet=0, const int* conn=0, bool addNode=false) : _nbSplits(nbTet), _nbCorners(4), _connectivity(conn), _baryNode(addNode), _ownConn(false) {} ~TSplitMethod() { if ( _ownConn ) delete [] _connectivity; _connectivity = 0; } + TSplitMethod(const TSplitMethod &splitMethod) + : _nbSplits(splitMethod._nbSplits), + _nbCorners(splitMethod._nbCorners), + _baryNode(splitMethod._baryNode), + _ownConn(splitMethod._ownConn), + _faceBaryNode(splitMethod._faceBaryNode) + { + _connectivity = splitMethod._connectivity; + const_cast(splitMethod)._connectivity = nullptr; + const_cast(splitMethod)._ownConn = false; + } bool hasFacet( const TTriangleFacet& facet ) const { if ( _nbCorners == 4 ) @@ -1904,7 +2272,7 @@ namespace } for ( int variant = 0; variant < nbVariants && method._nbSplits == 0; ++variant ) { - // check method compliancy with adjacent tetras, + // check method compliance with adjacent tetras, // all found splits must be among facets of tetras described by this method method = TSplitMethod( nbTet, connVariants[variant] ); if ( hasAdjacentSplits && method._nbSplits > 0 ) @@ -2035,6 +2403,8 @@ namespace const int methodFlags, const int facetToSplit) { + TSplitMethod method; + // order of facets in HEX according to SMDS_VolumeTool::Hexa_F : // B, T, L, B, R, F const int iF = ( facetToSplit < 2 ) ? 0 : 1 + ( facetToSplit-2 ) % 2; // [0,1,2] @@ -2065,12 +2435,12 @@ namespace to4methods[iF]._nbSplits = 4; to4methods[iF]._nbCorners = 6; } - return to4methods[iF]; + method = to4methods[iF]; + to4methods[iF]._connectivity = method._connectivity; // as copy ctor resets _connectivity + return method; } // else if ( methodFlags == HEXA_TO_2_PRISMS ) - TSplitMethod method; - const int iQ = vol.Element()->IsQuadratic() ? 2 : 1; const int nbVariants = 2, nbSplits = 2; @@ -2103,7 +2473,7 @@ namespace // there are adjacent prism for ( int variant = 0; variant < nbVariants; ++variant ) { - // check method compliancy with adjacent prisms, + // check method compliance with adjacent prisms, // the found prism facets must be among facets of prisms described by current method method._nbSplits = nbSplits; method._nbCorners = 6; @@ -2182,7 +2552,7 @@ namespace */ //======================================================================= - struct TVolumeFaceKey: pair< pair< int, int>, pair< int, int> > + struct TVolumeFaceKey: pair< pair< smIdType, smIdType>, pair< smIdType, smIdType> > { TVolumeFaceKey( SMDS_VolumeTool& vol, int iF ) { @@ -3096,10 +3466,10 @@ public: :myMesh( theMesh ), myMaxID( theMesh->MaxNodeID() + 1) {} - long GetLinkID (const SMDS_MeshNode * n1, + smIdType GetLinkID (const SMDS_MeshNode * n1, const SMDS_MeshNode * n2) const { - return ( Min(n1->GetID(),n2->GetID()) * myMaxID + Max(n1->GetID(),n2->GetID())); + return ( std::min(n1->GetID(),n2->GetID()) * myMaxID + std::max(n1->GetID(),n2->GetID())); } bool GetNodes (const long theLinkID, @@ -3482,6 +3852,68 @@ void SMESH_MeshEditor::GetLinkedNodes( const SMDS_MeshNode* theNode, } } +//======================================================================= +//function : averageBySurface +//purpose : Auxiliar function to treat properly nodes in periodic faces in the laplacian smoother +//======================================================================= +void averageBySurface( const Handle(Geom_Surface)& theSurface, const SMDS_MeshNode* refNode, + TIDSortedElemSet& nodeSet, map< const SMDS_MeshNode*, gp_XY* >& theUVMap, double * coord ) +{ + if ( theSurface.IsNull() ) + { + TIDSortedElemSet::iterator nodeSetIt = nodeSet.begin(); + for ( ; nodeSetIt != nodeSet.end(); nodeSetIt++ ) + { + const SMDS_MeshNode* node = cast2Node(*nodeSetIt); + coord[0] += node->X(); + coord[1] += node->Y(); + coord[2] += node->Z(); + } + } + else + { + Standard_Real Umin,Umax,Vmin,Vmax; + theSurface->Bounds( Umin, Umax, Vmin, Vmax ); + ASSERT( theUVMap.find( refNode ) != theUVMap.end() ); + gp_XY* nodeUV = theUVMap[ refNode ]; + Standard_Real uref = nodeUV->X(); + Standard_Real vref = nodeUV->Y(); + + TIDSortedElemSet::iterator nodeSetIt = nodeSet.begin(); + for ( ; nodeSetIt != nodeSet.end(); nodeSetIt++ ) + { + const SMDS_MeshNode* node = cast2Node(*nodeSetIt); + ASSERT( theUVMap.find( node ) != theUVMap.end() ); + gp_XY* uv = theUVMap[ node ]; + + if ( theSurface->IsUPeriodic() || theSurface->IsVPeriodic() ) + { + Standard_Real u = uv->X(); + Standard_Real v = uv->Y(); + Standard_Real uCorrected = u; + Standard_Real vCorrected = v; + bool isUTobeCorrected = (std::fabs( uref - u ) >= 0.7 * std::fabs( Umax - Umin )); + bool isVTobeCorrected = (std::fabs( vref - v ) >= 0.7 * std::fabs( Vmax - Vmin )); + + if( isUTobeCorrected ) + uCorrected = uref > u ? Umax + std::fabs(Umin - u) : Umin - std::fabs(Umax - u); + + if( isVTobeCorrected ) + vCorrected = vref > v ? Vmax + std::fabs(Vmin - v) : Vmin - std::fabs(Vmax - v); + + coord[0] += uCorrected; + coord[1] += vCorrected; + + } + else + { + coord[0] += uv->X(); + coord[1] += uv->Y(); + } + } + } +} + //======================================================================= //function : laplacianSmooth //purpose : pulls theNode toward the center of surrounding nodes directly @@ -3498,26 +3930,14 @@ void laplacianSmooth(const SMDS_MeshNode* theNode, SMESH_MeshEditor::GetLinkedNodes( theNode, nodeSet, SMDSAbs_Face ); // compute new coodrs + double coord[] = { 0., 0., 0. }; + + averageBySurface( theSurface, theNode, nodeSet, theUVMap, coord ); - double coord[] = { 0., 0., 0. }; - TIDSortedElemSet::iterator nodeSetIt = nodeSet.begin(); - for ( ; nodeSetIt != nodeSet.end(); nodeSetIt++ ) { - const SMDS_MeshNode* node = cast2Node(*nodeSetIt); - if ( theSurface.IsNull() ) { // smooth in 3D - coord[0] += node->X(); - coord[1] += node->Y(); - coord[2] += node->Z(); - } - else { // smooth in 2D - ASSERT( theUVMap.find( node ) != theUVMap.end() ); - gp_XY* uv = theUVMap[ node ]; - coord[0] += uv->X(); - coord[1] += uv->Y(); - } - } int nbNodes = nodeSet.size(); if ( !nbNodes ) return; + coord[0] /= nbNodes; coord[1] /= nbNodes; @@ -3527,7 +3947,7 @@ void laplacianSmooth(const SMDS_MeshNode* theNode, gp_Pnt p3d = theSurface->Value( coord[0], coord[1] ); coord[0] = p3d.X(); coord[1] = p3d.Y(); - coord[2] = p3d.Z(); + coord[2] = p3d.Z(); } else coord[2] /= nbNodes; @@ -3537,6 +3957,72 @@ void laplacianSmooth(const SMDS_MeshNode* theNode, const_cast< SMDS_MeshNode* >( theNode )->setXYZ(coord[0],coord[1],coord[2]); } +//======================================================================= +//function : correctTheValue +//purpose : Given a boundaries of parametric space determine if the node coordinate (u,v) need correction +// based on the reference coordinate (uref,vref) +//======================================================================= +void correctTheValue( Standard_Real Umax, Standard_Real Umin, Standard_Real Vmax, Standard_Real Vmin, + Standard_Real uref, Standard_Real vref, Standard_Real &u, Standard_Real &v ) +{ + bool isUTobeCorrected = (std::fabs( uref - u ) >= 0.7 * std::fabs( Umax - Umin )); + bool isVTobeCorrected = (std::fabs( vref - v ) >= 0.7 * std::fabs( Vmax - Vmin )); + if ( isUTobeCorrected ) + u = std::fabs(u-Umin) < 1e-7 ? Umax : Umin; + if ( isVTobeCorrected ) + v = std::fabs(v-Vmin) < 1e-7 ? Vmax : Vmin; +} + +//======================================================================= +//function : averageByElement +//purpose : Auxiliar function to treat properly nodes in periodic faces in the centroidal smoother +//======================================================================= +void averageByElement( const Handle(Geom_Surface)& theSurface, const SMDS_MeshNode* refNode, const SMDS_MeshElement* elem, + map< const SMDS_MeshNode*, gp_XY* >& theUVMap, SMESH::Controls::TSequenceOfXYZ& aNodePoints, + gp_XYZ& elemCenter ) +{ + int nn = elem->NbNodes(); + if(elem->IsQuadratic()) nn = nn/2; + int i=0; + SMDS_ElemIteratorPtr itN = elem->nodesIterator(); + Standard_Real Umin,Umax,Vmin,Vmax; + while ( i( itN->next() ); + i++; + gp_XYZ aP( aNode->X(), aNode->Y(), aNode->Z() ); + aNodePoints.push_back( aP ); + if ( !theSurface.IsNull() ) // smooth in 2D + { + ASSERT( theUVMap.find( aNode ) != theUVMap.end() ); + gp_XY* uv = theUVMap[ aNode ]; + + if ( theSurface->IsUPeriodic() || theSurface->IsVPeriodic() ) + { + theSurface->Bounds( Umin, Umax, Vmin, Vmax ); + Standard_Real u = uv->X(); + Standard_Real v = uv->Y(); + bool isSingularPoint = std::fabs(u - Umin) < 1e-7 || std::fabs(v - Vmin) < 1e-7 || std::fabs(u - Umax) < 1e-7 || std::fabs( v - Vmax ) < 1e-7; + if ( !isSingularPoint ) + { + aP.SetCoord( uv->X(), uv->Y(), 0. ); + } + else + { + gp_XY* refPoint = theUVMap[ refNode ]; + Standard_Real uref = refPoint->X(); + Standard_Real vref = refPoint->Y(); + correctTheValue( Umax, Umin, Vmax, Vmin, uref, vref, u, v ); + aP.SetCoord( u, v, 0. ); + } + } + else + aP.SetCoord( uv->X(), uv->Y(), 0. ); + } + elemCenter += aP; + } +} + //======================================================================= //function : centroidalSmooth //purpose : pulls theNode toward the element-area-weighted centroid of the @@ -3551,48 +4037,46 @@ void centroidalSmooth(const SMDS_MeshNode* theNode, SMESH::Controls::Area anAreaFunc; double totalArea = 0.; int nbElems = 0; - // compute new XYZ - + bool notToMoveNode = false; + // Do not correct singular nodes + if ( !theSurface.IsNull() && (theSurface->IsUPeriodic() || theSurface->IsVPeriodic()) ) + { + Standard_Real Umin,Umax,Vmin,Vmax; + theSurface->Bounds( Umin, Umax, Vmin, Vmax ); + gp_XY* uv = theUVMap[ theNode ]; + Standard_Real u = uv->X(); + Standard_Real v = uv->Y(); + notToMoveNode = std::fabs(u - Umin) < 1e-7 || std::fabs(v - Vmin) < 1e-7 || std::fabs(u - Umax) < 1e-7 || std::fabs( v - Vmax ) < 1e-7; + } + SMDS_ElemIteratorPtr elemIt = theNode->GetInverseElementIterator(SMDSAbs_Face); - while ( elemIt->more() ) + while ( elemIt->more() && !notToMoveNode ) { const SMDS_MeshElement* elem = elemIt->next(); nbElems++; gp_XYZ elemCenter(0.,0.,0.); SMESH::Controls::TSequenceOfXYZ aNodePoints; - SMDS_ElemIteratorPtr itN = elem->nodesIterator(); int nn = elem->NbNodes(); if(elem->IsQuadratic()) nn = nn/2; - int i=0; - //while ( itN->more() ) { - while ( i( itN->next() ); - i++; - gp_XYZ aP( aNode->X(), aNode->Y(), aNode->Z() ); - aNodePoints.push_back( aP ); - if ( !theSurface.IsNull() ) { // smooth in 2D - ASSERT( theUVMap.find( aNode ) != theUVMap.end() ); - gp_XY* uv = theUVMap[ aNode ]; - aP.SetCoord( uv->X(), uv->Y(), 0. ); - } - elemCenter += aP; - } + averageByElement( theSurface, theNode, elem, theUVMap, aNodePoints, elemCenter ); + double elemArea = anAreaFunc.GetValue( aNodePoints ); totalArea += elemArea; elemCenter /= nn; aNewXYZ += elemCenter * elemArea; } aNewXYZ /= totalArea; - if ( !theSurface.IsNull() ) { + + if ( !theSurface.IsNull() && !notToMoveNode ) { theUVMap[ theNode ]->SetCoord( aNewXYZ.X(), aNewXYZ.Y() ); aNewXYZ = theSurface->Value( aNewXYZ.X(), aNewXYZ.Y() ).XYZ(); } // move node - - const_cast< SMDS_MeshNode* >( theNode )->setXYZ(aNewXYZ.X(),aNewXYZ.Y(),aNewXYZ.Z()); + if ( !notToMoveNode ) + const_cast< SMDS_MeshNode* >( theNode )->setXYZ(aNewXYZ.X(),aNewXYZ.Y(),aNewXYZ.Z()); } //======================================================================= @@ -3606,7 +4090,7 @@ static bool getClosestUV (Extrema_GenExtPS& projector, { projector.Perform( point ); if ( projector.IsDone() ) { - double u, v, minVal = DBL_MAX; + double u = 0, v = 0, minVal = DBL_MAX; for ( int i = projector.NbExt(); i > 0; i-- ) if ( projector.SquareDistance( i ) < minVal ) { minVal = projector.SquareDistance( i ); @@ -3895,7 +4379,8 @@ void SMESH_MeshEditor::Smooth (TIDSortedElemSet & theElems, if ( !BRep_Tool::IsClosed( edge, face )) continue; SMESHDS_SubMesh* sm = aMesh->MeshElements( edge ); - if ( !sm ) continue; + if ( !sm ) + continue; // find out which parameter varies for a node on seam double f,l; gp_Pnt2d uv1, uv2; @@ -4269,7 +4754,7 @@ void SMESH_MeshEditor::sweepElement(const SMDS_MeshElement* elem, std::swap( itNN[0], itNN[1] ); std::swap( prevNod[0], prevNod[1] ); std::swap( nextNod[0], nextNod[1] ); - std::swap( isSingleNode[0], isSingleNode[1] ); + std::vector::swap(isSingleNode[0], isSingleNode[1]); if ( nbSame > 0 ) sames[0] = 1 - sames[0]; iNotSameNode = 1 - iNotSameNode; @@ -4852,13 +5337,12 @@ void SMESH_MeshEditor::makeWalls (TNodeOfNodeListMap & mapNewNodes, srcEdges.push_back(aMesh->FindEdge (commonNodes[0],commonNodes[1],commonNodes[2])); else srcEdges.push_back(aMesh->FindEdge (commonNodes[0],commonNodes[1])); -#ifdef _DEBUG_ - if ( !srcEdges.back() ) + + if (SALOME::VerbosityActivated() && !srcEdges.back()) { cout << "SMESH_MeshEditor::makeWalls(), no source edge found for a free face #" - << iF << " of volume #" << vTool.ID() << endl; + << iF << " of volume #" << vTool.ID() << endl; } -#endif } } if ( freeInd.empty() ) @@ -5184,6 +5668,7 @@ SMESH_MeshEditor::RotationSweep(TIDSortedElemSet theElemSets[2], SMESH_MeshEditor::ExtrusParam::ExtrusParam( const gp_Vec& theStep, const int theNbSteps, const std::list& theScales, + const std::list& theAngles, const gp_XYZ* theBasePoint, const int theFlags, const double theTolerance): @@ -5198,32 +5683,52 @@ SMESH_MeshEditor::ExtrusParam::ExtrusParam( const gp_Vec& theStep, for (int i=1; i<=theNbSteps; i++ ) mySteps->Append( stepSize ); - int nbScales = theScales.size(); - if ( nbScales > 0 ) + if ( !theScales.empty() ) { - if ( IsLinearVariation() && nbScales < theNbSteps ) + if ( IsScaleVariation() && (int)theScales.size() < theNbSteps ) + linearScaleVariation( theNbSteps, const_cast< std::list& >( theScales )); + + // add medium scales + std::list::const_iterator s2 = theScales.begin(), s1 = s2++; + myScales.reserve( theNbSteps * 2 ); + myScales.push_back( 0.5 * ( *s1 + 1. )); + myScales.push_back( *s1 ); + for ( ; s2 != theScales.end(); s1 = s2++ ) { - myScales.reserve( theNbSteps ); - std::list::const_iterator scale = theScales.begin(); - double prevScale = 1.0; - for ( int iSc = 1; scale != theScales.end(); ++scale, ++iSc ) - { - int iStep = int( iSc / double( nbScales ) * theNbSteps + 0.5 ); - int stDelta = Max( 1, iStep - myScales.size()); - double scDelta = ( *scale - prevScale ) / stDelta; - for ( int iStep = 0; iStep < stDelta; ++iStep ) - { - myScales.push_back( prevScale + scDelta ); - prevScale = myScales.back(); - } - prevScale = *scale; - } + myScales.push_back( 0.5 * ( *s1 + *s2 )); + myScales.push_back( *s2 ); } - else + } + + if ( !theAngles.empty() ) + { + std::list& angles = const_cast< std::list& >( theAngles ); + if ( IsAngleVariation() && (int)theAngles.size() < theNbSteps ) + linearAngleVariation( theNbSteps, angles ); + + // accumulate angles + double angle = 0; + int nbAngles = 0; + std::list::iterator a1 = angles.begin(), a2; + for ( ; a1 != angles.end(); ++a1, ++nbAngles ) + { + angle += *a1; + *a1 = angle; + } + while ( nbAngles++ < theNbSteps ) + angles.push_back( angles.back() ); + + // add medium angles + a2 = angles.begin(), a1 = a2++; + myAngles.push_back( 0.5 * *a1 ); + myAngles.push_back( *a1 ); + for ( ; a2 != angles.end(); a1 = a2++ ) { - myScales.assign( theScales.begin(), theScales.end() ); + myAngles.push_back( 0.5 * ( *a1 + *a2 )); + myAngles.push_back( *a2 ); } } + if ( theBasePoint ) { myBaseP = *theBasePoint; @@ -5294,6 +5799,41 @@ SMESH_MeshEditor::ExtrusParam::ExtrusParam( const double theStepSize, } } +//======================================================================= +//function : ExtrusParam +//purpose : for extrusion along path +//======================================================================= + +SMESH_MeshEditor::ExtrusParam::ExtrusParam( const std::vector< PathPoint >& thePoints, + const gp_Pnt* theBasePoint, + const std::list& theScales, + const bool theMakeGroups ) + : myBaseP( Precision::Infinite(), 0, 0 ), + myFlags( EXTRUSION_FLAG_BOUNDARY | ( theMakeGroups ? EXTRUSION_FLAG_GROUPS : 0 )), + myPathPoints( thePoints ) +{ + if ( theBasePoint ) + { + myBaseP = theBasePoint->XYZ(); + } + + if ( !theScales.empty() ) + { + // add medium scales + std::list::const_iterator s2 = theScales.begin(), s1 = s2++; + myScales.reserve( thePoints.size() * 2 ); + myScales.push_back( 0.5 * ( 1. + *s1 )); + myScales.push_back( *s1 ); + for ( ; s2 != theScales.end(); s1 = s2++ ) + { + myScales.push_back( 0.5 * ( *s1 + *s2 )); + myScales.push_back( *s2 ); + } + } + + myMakeNodesFun = & SMESH_MeshEditor::ExtrusParam::makeNodesAlongTrack; +} + //======================================================================= //function : ExtrusParam::SetElementsToUse //purpose : stores elements to use for extrusion by normal, depending on @@ -5399,38 +5939,37 @@ makeNodesByDir( SMESHDS_Mesh* mesh, newNodes.push_back( newNode ); } - if ( !myScales.empty() ) + if ( !myScales.empty() || !myAngles.empty() ) { - if ( makeMediumNodes && myMediumScales.empty() ) - { - myMediumScales.resize( myScales.size() ); - double prevFactor = 1.; - for ( size_t i = 0; i < myScales.size(); ++i ) - { - myMediumScales[i] = 0.5 * ( prevFactor + myScales[i] ); - prevFactor = myScales[i]; - } - } - typedef std::vector::iterator ScaleIt; - ScaleIt scales[] = { myScales.begin(), myMediumScales.begin() }; - - size_t iSc = 0, nbScales = myScales.size() + myMediumScales.size(); + gp_XYZ center = myBaseP; + gp_Ax1 ratationAxis( center, myDir ); + gp_Trsf rotation; - gp_XYZ center = myBaseP; std::list::iterator nIt = newNodes.begin(); - size_t iN = 0; - for ( beginStepIter( makeMediumNodes ); moreSteps() && ( iN < nbScales ); ++nIt, ++iN ) + size_t i = !makeMediumNodes; + for ( beginStepIter( makeMediumNodes ); + moreSteps(); + ++nIt, i += 1 + !makeMediumNodes ) { center += myDir.XYZ() * nextStep(); - iSc += int( makeMediumNodes ); - ScaleIt& scale = scales[ iSc % 2 ]; - gp_XYZ xyz = SMESH_NodeXYZ( *nIt ); - xyz = ( *scale * ( xyz - center )) + center; - mesh->MoveNode( *nIt, xyz.X(), xyz.Y(), xyz.Z() ); - - ++scale; + bool moved = false; + if ( i < myScales.size() ) + { + xyz = ( myScales[i] * ( xyz - center )) + center; + moved = true; + } + if ( !myAngles.empty() ) + { + rotation.SetRotation( ratationAxis, myAngles[i] ); + rotation.Transforms( xyz ); + moved = true; + } + if ( moved ) + mesh->MoveNode( *nIt, xyz.X(), xyz.Y(), xyz.Z() ); + else + break; } } return nbNodes; @@ -5588,13 +6127,107 @@ makeNodesByNormal2D( SMESHDS_Mesh* mesh, //======================================================================= int SMESH_MeshEditor::ExtrusParam:: -makeNodesByNormal1D( SMESHDS_Mesh* mesh, +makeNodesByNormal1D( SMESHDS_Mesh* /*mesh*/, + const SMDS_MeshNode* /*srcNode*/, + std::list & /*newNodes*/, + const bool /*makeMediumNodes*/) +{ + throw SALOME_Exception("Extrusion 1D by Normal not implemented"); + return 0; +} + +//======================================================================= +//function : ExtrusParam::makeNodesAlongTrack +//purpose : create nodes for extrusion along path +//======================================================================= + +int SMESH_MeshEditor::ExtrusParam:: +makeNodesAlongTrack( SMESHDS_Mesh* mesh, const SMDS_MeshNode* srcNode, std::list & newNodes, const bool makeMediumNodes) { - throw SALOME_Exception("Extrusion 1D by Normal not implemented"); - return 0; + const Standard_Real aTolAng=1.e-4; + + gp_Pnt aV0x = myBaseP; + gp_Pnt aPN0 = SMESH_NodeXYZ( srcNode ); + + const PathPoint& aPP0 = myPathPoints[0]; + gp_Pnt aP0x = aPP0.myPnt; + gp_Dir aDT0x= aPP0.myTgt; + + std::vector< gp_Pnt > centers; + centers.reserve( NbSteps() * 2 ); + + gp_Trsf aTrsf, aTrsfRot, aTrsfRotT1T0; + + for ( size_t j = 1; j < myPathPoints.size(); ++j ) + { + const PathPoint& aPP = myPathPoints[j]; + const gp_Pnt& aP1x = aPP.myPnt; + const gp_Dir& aDT1x = aPP.myTgt; + + // Translation + gp_Vec aV01x( aP0x, aP1x ); + aTrsf.SetTranslation( aV01x ); + gp_Pnt aV1x = aV0x.Transformed( aTrsf ); + gp_Pnt aPN1 = aPN0.Transformed( aTrsf ); + + // rotation 1 [ T1,T0 ] + Standard_Real aAngleT1T0 = -aDT1x.Angle( aDT0x ); + if ( fabs( aAngleT1T0 ) > aTolAng ) + { + gp_Dir aDT1T0 = aDT1x ^ aDT0x; + aTrsfRotT1T0.SetRotation( gp_Ax1( aV1x, aDT1T0 ), aAngleT1T0 ); + + aPN1 = aPN1.Transformed( aTrsfRotT1T0 ); + } + + // rotation 2 + if ( aPP.myAngle != 0. ) + { + aTrsfRot.SetRotation( gp_Ax1( aV1x, aDT1x ), aPP.myAngle ); + aPN1 = aPN1.Transformed( aTrsfRot ); + } + + // make new node + if ( makeMediumNodes ) + { + // create additional node + gp_XYZ midP = 0.5 * ( aPN1.XYZ() + aPN0.XYZ() ); + const SMDS_MeshNode* newNode = mesh->AddNode( midP.X(), midP.Y(), midP.Z() ); + newNodes.push_back( newNode ); + + } + const SMDS_MeshNode* newNode = mesh->AddNode( aPN1.X(), aPN1.Y(), aPN1.Z() ); + newNodes.push_back( newNode ); + + centers.push_back( 0.5 * ( aV0x.XYZ() + aV1x.XYZ() )); + centers.push_back( aV1x ); + + aPN0 = aPN1; + aP0x = aP1x; + aV0x = aV1x; + aDT0x = aDT1x; + } + + // scale + if ( !myScales.empty() ) + { + gp_Trsf aTrsfScale; + std::list::iterator node = newNodes.begin(); + for ( size_t i = !makeMediumNodes; + i < myScales.size() && node != newNodes.end(); + i += ( 1 + !makeMediumNodes ), ++node ) + { + aTrsfScale.SetScale( centers[ i ], myScales[ i ] ); + gp_Pnt aN = SMESH_NodeXYZ( *node ); + gp_Pnt aP = aN.Transformed( aTrsfScale ); + mesh->MoveNode( *node, aP.X(), aP.Y(), aP.Z() ); + } + } + + return myPathPoints.size() + makeMediumNodes * ( myPathPoints.size() - 2 ); } //======================================================================= @@ -5610,10 +6243,31 @@ SMESH_MeshEditor::ExtrusionSweep (TIDSortedElemSet theElems[2], const int theFlags, const double theTolerance) { - ExtrusParam aParams( theStep, theNbSteps, std::list(), 0, theFlags, theTolerance ); + std::list dummy; + ExtrusParam aParams( theStep, theNbSteps, dummy, dummy, 0, + theFlags, theTolerance ); return ExtrusionSweep( theElems, aParams, newElemsMap ); } +namespace +{ + +//======================================================================= +//function : getOriFactor +//purpose : Return -1 or 1 depending on if order of given nodes corresponds to +// edge curve orientation +//======================================================================= + + double getOriFactor( const TopoDS_Edge& edge, + const SMDS_MeshNode* n1, + const SMDS_MeshNode* n2, + SMESH_MesherHelper& helper) + { + double u1 = helper.GetNodeU( edge, n1, n2 ); + double u2 = helper.GetNodeU( edge, n2, n1 ); + return u1 < u2 ? 1. : -1.; + } +} //======================================================================= //function : ExtrusionSweep @@ -5662,11 +6316,11 @@ SMESH_MeshEditor::ExtrusionSweep (TIDSortedElemSet theElemSets[2], newNodesItVec.reserve( nbNodes ); // loop on elem nodes - SMDS_ElemIteratorPtr itN = elem->nodesIterator(); + SMDS_NodeIteratorPtr itN = elem->nodeIterator(); while ( itN->more() ) { // check if a node has been already sweeped - const SMDS_MeshNode* node = cast2Node( itN->next() ); + const SMDS_MeshNode* node = itN->next(); TNodeOfNodeListMap::iterator nIt = mapNewNodes.insert( make_pair( node, list() )).first; list& listNewNodes = nIt->second; @@ -5700,6 +6354,11 @@ SMESH_MeshEditor::ExtrusionSweep (TIDSortedElemSet theElemSets[2], } else { + if ( theParams.ToMakeBoundary() ) + { + GetMeshDS()->Modified(); + throw SALOME_Exception( SMESH_Comment("Can't extrude node #") << node->GetID() ); + } break; // newNodesItVec will be shorter than nbNodes } } @@ -5727,726 +6386,184 @@ SMESH_MeshEditor::ExtrusionSweep (TIDSortedElemSet theElemSets[2], //======================================================================= SMESH_MeshEditor::Extrusion_Error SMESH_MeshEditor::ExtrusionAlongTrack (TIDSortedElemSet theElements[2], - SMESH_subMesh* theTrack, + SMESH_Mesh* theTrackMesh, + SMDS_ElemIteratorPtr theTrackIterator, const SMDS_MeshNode* theN1, - const bool theHasAngles, - list& theAngles, - const bool theLinearVariation, - const bool theHasRefPoint, - const gp_Pnt& theRefPoint, + std::list& theAngles, + const bool theAngleVariation, + std::list& theScales, + const bool theScaleVariation, + const gp_Pnt* theRefPoint, const bool theMakeGroups) { ClearLastCreated(); - int aNbE; - std::list aPrms; - TIDSortedElemSet::iterator itElem; - - gp_XYZ aGC; - TopoDS_Edge aTrackEdge; - TopoDS_Vertex aV1, aV2; - - SMDS_ElemIteratorPtr aItE; - SMDS_NodeIteratorPtr aItN; - SMDSAbs_ElementType aTypeE; - - TNodeOfNodeListMap mapNewNodes; - // 1. Check data - aNbE = theElements[0].size() + theElements[1].size(); - // nothing to do - if ( !aNbE ) + if ( theElements[0].empty() && theElements[1].empty() ) return EXTR_NO_ELEMENTS; - // 1.1 Track Pattern - ASSERT( theTrack ); - - SMESHDS_SubMesh* pSubMeshDS = theTrack->GetSubMeshDS(); - if ( !pSubMeshDS ) - return ExtrusionAlongTrack( theElements, theTrack->GetFather(), theN1, - theHasAngles, theAngles, theLinearVariation, - theHasRefPoint, theRefPoint, theMakeGroups ); - - aItE = pSubMeshDS->GetElements(); - while ( aItE->more() ) { - const SMDS_MeshElement* pE = aItE->next(); - aTypeE = pE->GetType(); - // Pattern must contain links only - if ( aTypeE != SMDSAbs_Edge ) - return EXTR_PATH_NOT_EDGE; - } - - list fullList; - - const TopoDS_Shape& aS = theTrack->GetSubShape(); - // Sub-shape for the Pattern must be an Edge or Wire - if( aS.ShapeType() == TopAbs_EDGE ) { - aTrackEdge = TopoDS::Edge( aS ); - // the Edge must not be degenerated - if ( SMESH_Algo::isDegenerated( aTrackEdge ) ) - return EXTR_BAD_PATH_SHAPE; - TopExp::Vertices( aTrackEdge, aV1, aV2 ); - aItN = theTrack->GetFather()->GetSubMesh( aV1 )->GetSubMeshDS()->GetNodes(); - const SMDS_MeshNode* aN1 = aItN->next(); - aItN = theTrack->GetFather()->GetSubMesh( aV2 )->GetSubMeshDS()->GetNodes(); - const SMDS_MeshNode* aN2 = aItN->next(); - // starting node must be aN1 or aN2 - if ( !( aN1 == theN1 || aN2 == theN1 ) ) - return EXTR_BAD_STARTING_NODE; - aItN = pSubMeshDS->GetNodes(); - while ( aItN->more() ) { - const SMDS_MeshNode* pNode = aItN->next(); - SMDS_EdgePositionPtr pEPos = pNode->GetPosition(); - double aT = pEPos->GetUParameter(); - aPrms.push_back( aT ); - } - //Extrusion_Error err = - makeEdgePathPoints(aPrms, aTrackEdge, (aN1==theN1), fullList); - } else if( aS.ShapeType() == TopAbs_WIRE ) { - list< SMESH_subMesh* > LSM; - TopTools_SequenceOfShape Edges; - SMESH_subMeshIteratorPtr itSM = theTrack->getDependsOnIterator(false,true); - while(itSM->more()) { - SMESH_subMesh* SM = itSM->next(); - LSM.push_back(SM); - const TopoDS_Shape& aS = SM->GetSubShape(); - Edges.Append(aS); - } - list< list > LLPPs; - int startNid = theN1->GetID(); - TColStd_MapOfInteger UsedNums; - - int NbEdges = Edges.Length(); - int i = 1; - for(; i<=NbEdges; i++) { - int k = 0; - list< SMESH_subMesh* >::iterator itLSM = LSM.begin(); - for(; itLSM!=LSM.end(); itLSM++) { - k++; - if(UsedNums.Contains(k)) continue; - aTrackEdge = TopoDS::Edge( Edges.Value(k) ); - SMESH_subMesh* locTrack = *itLSM; - SMESHDS_SubMesh* locMeshDS = locTrack->GetSubMeshDS(); - TopExp::Vertices( aTrackEdge, aV1, aV2 ); - aItN = locTrack->GetFather()->GetSubMesh(aV1)->GetSubMeshDS()->GetNodes(); - const SMDS_MeshNode* aN1 = aItN->next(); - aItN = locTrack->GetFather()->GetSubMesh(aV2)->GetSubMeshDS()->GetNodes(); - const SMDS_MeshNode* aN2 = aItN->next(); - // starting node must be aN1 or aN2 - if ( !( aN1->GetID() == startNid || aN2->GetID() == startNid ) ) continue; - // 2. Collect parameters on the track edge - aPrms.clear(); - aItN = locMeshDS->GetNodes(); - while ( aItN->more() ) { - const SMDS_MeshNode* pNode = aItN->next(); - SMDS_EdgePositionPtr pEPos = pNode->GetPosition(); - double aT = pEPos->GetUParameter(); - aPrms.push_back( aT ); - } - list LPP; - //Extrusion_Error err = - makeEdgePathPoints(aPrms, aTrackEdge,(aN1->GetID()==startNid), LPP); - LLPPs.push_back(LPP); - UsedNums.Add(k); - // update startN for search following edge - if( aN1->GetID() == startNid ) startNid = aN2->GetID(); - else startNid = aN1->GetID(); - break; - } - } - list< list >::iterator itLLPP = LLPPs.begin(); - list firstList = *itLLPP; - list::iterator itPP = firstList.begin(); - for(; itPP!=firstList.end(); itPP++) { - fullList.push_back( *itPP ); - } - SMESH_MeshEditor_PathPoint PP1 = fullList.back(); - fullList.pop_back(); - itLLPP++; - for(; itLLPP!=LLPPs.end(); itLLPP++) { - list currList = *itLLPP; - itPP = currList.begin(); - SMESH_MeshEditor_PathPoint PP2 = currList.front(); - gp_Dir D1 = PP1.Tangent(); - gp_Dir D2 = PP2.Tangent(); - gp_Dir Dnew( gp_Vec( (D1.X()+D2.X())/2, (D1.Y()+D2.Y())/2, - (D1.Z()+D2.Z())/2 ) ); - PP1.SetTangent(Dnew); - fullList.push_back(PP1); - itPP++; - for(; itPP!=firstList.end(); itPP++) { - fullList.push_back( *itPP ); - } - PP1 = fullList.back(); - fullList.pop_back(); - } - // if wire not closed - fullList.push_back(PP1); - // else ??? - } - else { - return EXTR_BAD_PATH_SHAPE; - } - - return makeExtrElements(theElements, fullList, theHasAngles, theAngles, theLinearVariation, - theHasRefPoint, theRefPoint, theMakeGroups); -} - - -//======================================================================= -//function : ExtrusionAlongTrack -//purpose : -//======================================================================= -SMESH_MeshEditor::Extrusion_Error -SMESH_MeshEditor::ExtrusionAlongTrack (TIDSortedElemSet theElements[2], - SMESH_Mesh* theTrack, - const SMDS_MeshNode* theN1, - const bool theHasAngles, - list& theAngles, - const bool theLinearVariation, - const bool theHasRefPoint, - const gp_Pnt& theRefPoint, - const bool theMakeGroups) -{ - ClearLastCreated(); - - int aNbE; - std::list aPrms; - TIDSortedElemSet::iterator itElem; - - gp_XYZ aGC; - TopoDS_Edge aTrackEdge; - TopoDS_Vertex aV1, aV2; - - SMDS_ElemIteratorPtr aItE; - SMDS_NodeIteratorPtr aItN; - SMDSAbs_ElementType aTypeE; + ASSERT( theTrackMesh ); + if ( ! theTrackIterator || !theTrackIterator->more() ) + return EXTR_NO_ELEMENTS; - TNodeOfNodeListMap mapNewNodes; + // 2. Get ordered nodes + SMESH_MeshAlgos::TElemGroupVector branchEdges; + SMESH_MeshAlgos::TNodeGroupVector branchNods; + SMESH_MeshAlgos::Get1DBranches( theTrackIterator, branchEdges, branchNods, theN1 ); + if ( branchEdges.empty() ) + return EXTR_PATH_NOT_EDGE; - // 1. Check data - aNbE = theElements[0].size() + theElements[1].size(); - // nothing to do - if ( !aNbE ) - return EXTR_NO_ELEMENTS; + if ( branchEdges.size() > 1 ) + return EXTR_BAD_PATH_SHAPE; - // 1.1 Track Pattern - ASSERT( theTrack ); - - SMESHDS_Mesh* pMeshDS = theTrack->GetMeshDS(); - - aItE = pMeshDS->elementsIterator(); - while ( aItE->more() ) { - const SMDS_MeshElement* pE = aItE->next(); - aTypeE = pE->GetType(); - // Pattern must contain links only - if ( aTypeE != SMDSAbs_Edge ) - return EXTR_PATH_NOT_EDGE; - } - - list fullList; - - const TopoDS_Shape& aS = theTrack->GetShapeToMesh(); - - if ( !theTrack->HasShapeToMesh() ) { - //Mesh without shape - const SMDS_MeshNode* currentNode = NULL; - const SMDS_MeshNode* prevNode = theN1; - std::vector aNodesList; - aNodesList.push_back(theN1); - int nbEdges = 0, conn=0; - const SMDS_MeshElement* prevElem = NULL; - const SMDS_MeshElement* currentElem = NULL; - int totalNbEdges = theTrack->NbEdges(); - SMDS_ElemIteratorPtr nIt; - - //check start node - if( !theTrack->GetMeshDS()->Contains( theN1 )) { - return EXTR_BAD_STARTING_NODE; - } - - conn = nbEdgeConnectivity(theN1); - if( conn != 1 ) - return EXTR_PATH_NOT_EDGE; - - aItE = theN1->GetInverseElementIterator(); - prevElem = aItE->next(); - currentElem = prevElem; - //Get all nodes - if(totalNbEdges == 1 ) { - nIt = currentElem->nodesIterator(); - currentNode = static_cast(nIt->next()); - if(currentNode == prevNode) - currentNode = static_cast(nIt->next()); - aNodesList.push_back(currentNode); - } else { - nIt = currentElem->nodesIterator(); - while( nIt->more() ) { - currentNode = static_cast(nIt->next()); - if(currentNode == prevNode) - currentNode = static_cast(nIt->next()); - aNodesList.push_back(currentNode); - - //case of the closed mesh - if(currentNode == theN1) { - nbEdges++; - break; - } + std::vector< const SMDS_MeshNode* >& pathNodes = branchNods[0]; + std::vector< const SMDS_MeshElement* >& pathEdges = branchEdges[0]; + if ( pathNodes[0] != theN1 && pathNodes[1] != theN1 ) + return EXTR_BAD_STARTING_NODE; - conn = nbEdgeConnectivity(currentNode); - if(conn > 2) { - return EXTR_PATH_NOT_EDGE; - }else if( conn == 1 && nbEdges > 0 ) { - //End of the path - nbEdges++; - break; - }else { - prevNode = currentNode; - aItE = currentNode->GetInverseElementIterator(); - currentElem = aItE->next(); - if( currentElem == prevElem) - currentElem = aItE->next(); - nIt = currentElem->nodesIterator(); - prevElem = currentElem; - nbEdges++; - } - } - } - - if(nbEdges != totalNbEdges) - return EXTR_PATH_NOT_EDGE; - - TopTools_SequenceOfShape Edges; - list< list > LLPPs; - int startNid = theN1->GetID(); - for ( size_t i = 1; i < aNodesList.size(); i++ ) - { - gp_Pnt p1 = SMESH_NodeXYZ( aNodesList[i-1] ); - gp_Pnt p2 = SMESH_NodeXYZ( aNodesList[i] ); - TopoDS_Edge e = BRepBuilderAPI_MakeEdge( p1, p2 ); - list LPP; - aPrms.clear(); - makeEdgePathPoints(aPrms, e, (aNodesList[i-1]->GetID()==startNid), LPP); - LLPPs.push_back(LPP); - if ( aNodesList[i-1]->GetID() == startNid ) startNid = aNodesList[i ]->GetID(); - else startNid = aNodesList[i-1]->GetID(); - } - - list< list >::iterator itLLPP = LLPPs.begin(); - list firstList = *itLLPP; - list::iterator itPP = firstList.begin(); - for(; itPP!=firstList.end(); itPP++) { - fullList.push_back( *itPP ); - } - - SMESH_MeshEditor_PathPoint PP1 = fullList.back(); - SMESH_MeshEditor_PathPoint PP2; - fullList.pop_back(); - itLLPP++; - for(; itLLPP!=LLPPs.end(); itLLPP++) { - list currList = *itLLPP; - itPP = currList.begin(); - PP2 = currList.front(); - gp_Dir D1 = PP1.Tangent(); - gp_Dir D2 = PP2.Tangent(); - gp_Dir Dnew( 0.5 * ( D1.XYZ() + D2.XYZ() )); - PP1.SetTangent(Dnew); - fullList.push_back(PP1); - itPP++; - for(; itPP!=currList.end(); itPP++) { - fullList.push_back( *itPP ); - } - PP1 = fullList.back(); - fullList.pop_back(); - } - fullList.push_back(PP1); - - } // Sub-shape for the Pattern must be an Edge or Wire - else if ( aS.ShapeType() == TopAbs_EDGE ) - { - aTrackEdge = TopoDS::Edge( aS ); - // the Edge must not be degenerated - if ( SMESH_Algo::isDegenerated( aTrackEdge ) ) - return EXTR_BAD_PATH_SHAPE; - TopExp::Vertices( aTrackEdge, aV1, aV2 ); - const SMDS_MeshNode* aN1 = SMESH_Algo::VertexNode( aV1, pMeshDS ); - const SMDS_MeshNode* aN2 = SMESH_Algo::VertexNode( aV2, pMeshDS ); - // starting node must be aN1 or aN2 - if ( !( aN1 == theN1 || aN2 == theN1 ) ) - return EXTR_BAD_STARTING_NODE; - aItN = pMeshDS->nodesIterator(); - while ( aItN->more() ) { - const SMDS_MeshNode* pNode = aItN->next(); - if( pNode==aN1 || pNode==aN2 ) continue; - SMDS_EdgePositionPtr pEPos = pNode->GetPosition(); - double aT = pEPos->GetUParameter(); - aPrms.push_back( aT ); - } - //Extrusion_Error err = - makeEdgePathPoints(aPrms, aTrackEdge, (aN1==theN1), fullList); - } - else if( aS.ShapeType() == TopAbs_WIRE ) { - list< SMESH_subMesh* > LSM; - TopTools_SequenceOfShape Edges; - TopExp_Explorer eExp(aS, TopAbs_EDGE); - for(; eExp.More(); eExp.Next()) { - TopoDS_Edge E = TopoDS::Edge( eExp.Current() ); - if( SMESH_Algo::isDegenerated(E) ) continue; - SMESH_subMesh* SM = theTrack->GetSubMesh(E); - if(SM) { - LSM.push_back(SM); - Edges.Append(E); - } - } - list< list > LLPPs; - TopoDS_Vertex aVprev; - TColStd_MapOfInteger UsedNums; - int NbEdges = Edges.Length(); - int i = 1; - for(; i<=NbEdges; i++) { - int k = 0; - list< SMESH_subMesh* >::iterator itLSM = LSM.begin(); - for(; itLSM!=LSM.end(); itLSM++) { - k++; - if(UsedNums.Contains(k)) continue; - aTrackEdge = TopoDS::Edge( Edges.Value(k) ); - SMESH_subMesh* locTrack = *itLSM; - SMESHDS_SubMesh* locMeshDS = locTrack->GetSubMeshDS(); - TopExp::Vertices( aTrackEdge, aV1, aV2 ); - bool aN1isOK = false, aN2isOK = false; - if ( aVprev.IsNull() ) { - // if previous vertex is not yet defined, it means that we in the beginning of wire - // and we have to find initial vertex corresponding to starting node theN1 - const SMDS_MeshNode* aN1 = SMESH_Algo::VertexNode( aV1, pMeshDS ); - const SMDS_MeshNode* aN2 = SMESH_Algo::VertexNode( aV2, pMeshDS ); - // starting node must be aN1 or aN2 - aN1isOK = ( aN1 && aN1 == theN1 ); - aN2isOK = ( aN2 && aN2 == theN1 ); - } - else { - // we have specified ending vertex of the previous edge on the previous iteration - // and we have just to check that it corresponds to any vertex in current segment - aN1isOK = aVprev.IsSame( aV1 ); - aN2isOK = aVprev.IsSame( aV2 ); - } - if ( !aN1isOK && !aN2isOK ) continue; - // 2. Collect parameters on the track edge - aPrms.clear(); - aItN = locMeshDS->GetNodes(); - while ( aItN->more() ) { - const SMDS_MeshNode* pNode = aItN->next(); - SMDS_EdgePositionPtr pEPos = pNode->GetPosition(); - double aT = pEPos->GetUParameter(); - aPrms.push_back( aT ); - } - list LPP; - //Extrusion_Error err = - makeEdgePathPoints(aPrms, aTrackEdge, aN1isOK, LPP); - LLPPs.push_back(LPP); - UsedNums.Add(k); - // update startN for search following edge - if ( aN1isOK ) aVprev = aV2; - else aVprev = aV1; - break; + if ( theTrackMesh->NbEdges( ORDER_QUADRATIC ) > 0 ) + { + // add medium nodes to pathNodes + std::vector< const SMDS_MeshNode* > pathNodes2; + std::vector< const SMDS_MeshElement* > pathEdges2; + pathNodes2.reserve( pathNodes.size() * 2 ); + pathEdges2.reserve( pathEdges.size() * 2 ); + for ( size_t i = 0; i < pathEdges.size(); ++i ) + { + pathNodes2.push_back( pathNodes[i] ); + pathEdges2.push_back( pathEdges[i] ); + if ( pathEdges[i]->IsQuadratic() ) + { + pathNodes2.push_back( pathEdges[i]->GetNode(2) ); + pathEdges2.push_back( pathEdges[i] ); } } - list< list >::iterator itLLPP = LLPPs.begin(); - list& firstList = *itLLPP; - fullList.splice( fullList.end(), firstList ); - - SMESH_MeshEditor_PathPoint PP1 = fullList.back(); - fullList.pop_back(); - itLLPP++; - for(; itLLPP!=LLPPs.end(); itLLPP++) { - list& currList = *itLLPP; - SMESH_MeshEditor_PathPoint PP2 = currList.front(); - gp_Dir D1 = PP1.Tangent(); - gp_Dir D2 = PP2.Tangent(); - gp_Dir Dnew( D1.XYZ() + D2.XYZ() ); - PP1.SetTangent(Dnew); - fullList.push_back(PP1); - fullList.splice( fullList.end(), currList, ++currList.begin(), currList.end() ); - PP1 = fullList.back(); - fullList.pop_back(); - } - // if wire not closed - fullList.push_back(PP1); - // else ??? - } - else { - return EXTR_BAD_PATH_SHAPE; + pathNodes2.push_back( pathNodes.back() ); + pathEdges.swap( pathEdges2 ); + pathNodes.swap( pathNodes2 ); } - return makeExtrElements(theElements, fullList, theHasAngles, theAngles, theLinearVariation, - theHasRefPoint, theRefPoint, theMakeGroups); -} + // 3. Get path data at pathNodes + std::vector< ExtrusParam::PathPoint > points( pathNodes.size() ); -//======================================================================= -//function : makeEdgePathPoints -//purpose : auxiliary for ExtrusionAlongTrack -//======================================================================= -SMESH_MeshEditor::Extrusion_Error -SMESH_MeshEditor::makeEdgePathPoints(std::list& aPrms, - const TopoDS_Edge& aTrackEdge, - bool FirstIsStart, - list& LPP) -{ - Standard_Real aTx1, aTx2, aL2, aTolVec, aTolVec2; - aTolVec=1.e-7; - aTolVec2=aTolVec*aTolVec; - double aT1, aT2; - TopoDS_Vertex aV1, aV2; - TopExp::Vertices( aTrackEdge, aV1, aV2 ); - aT1=BRep_Tool::Parameter( aV1, aTrackEdge ); - aT2=BRep_Tool::Parameter( aV2, aTrackEdge ); - // 2. Collect parameters on the track edge - aPrms.push_front( aT1 ); - aPrms.push_back( aT2 ); - // sort parameters - aPrms.sort(); - if( FirstIsStart ) { - if ( aT1 > aT2 ) { - aPrms.reverse(); - } - } - else { - if ( aT2 > aT1 ) { - aPrms.reverse(); - } - } - // 3. Path Points - SMESH_MeshEditor_PathPoint aPP; - Handle(Geom_Curve) aC3D = BRep_Tool::Curve( aTrackEdge, aTx1, aTx2 ); - std::list::iterator aItD = aPrms.begin(); - for(; aItD != aPrms.end(); ++aItD) { - double aT = *aItD; - gp_Pnt aP3D; - gp_Vec aVec; - aC3D->D1( aT, aP3D, aVec ); - aL2 = aVec.SquareMagnitude(); - if ( aL2 < aTolVec2 ) - return EXTR_CANT_GET_TANGENT; - gp_Dir aTgt( FirstIsStart ? aVec : -aVec ); - aPP.SetPnt( aP3D ); - aPP.SetTangent( aTgt ); - aPP.SetParameter( aT ); - LPP.push_back(aPP); - } - return EXTR_OK; -} + if ( theAngleVariation ) + linearAngleVariation( points.size()-1, theAngles ); + if ( theScaleVariation ) + linearScaleVariation( points.size()-1, theScales ); + theAngles.push_front( 0 ); // for the 1st point that is not transformed + std::list::iterator angle = theAngles.begin(); -//======================================================================= -//function : makeExtrElements -//purpose : auxiliary for ExtrusionAlongTrack -//======================================================================= -SMESH_MeshEditor::Extrusion_Error -SMESH_MeshEditor::makeExtrElements(TIDSortedElemSet theElemSets[2], - list& fullList, - const bool theHasAngles, - list& theAngles, - const bool theLinearVariation, - const bool theHasRefPoint, - const gp_Pnt& theRefPoint, - const bool theMakeGroups) -{ - const int aNbTP = fullList.size(); + SMESHDS_Mesh* pathMeshDS = theTrackMesh->GetMeshDS(); - // Angles - if( theHasAngles && !theAngles.empty() && theLinearVariation ) - linearAngleVariation(aNbTP-1, theAngles); + std::map< int, double > edgeID2OriFactor; // orientation of EDGEs + std::map< int, double >::iterator id2factor; + SMESH_MesherHelper pathHelper( *theTrackMesh ); + gp_Pnt p; gp_Vec tangent; + const double tol2 = gp::Resolution() * gp::Resolution(); - // fill vector of path points with angles - vector aPPs; - list::iterator itPP = fullList.begin(); - list::iterator itAngles = theAngles.begin(); - aPPs.push_back( *itPP++ ); - for( ; itPP != fullList.end(); itPP++) { - aPPs.push_back( *itPP ); - if ( theHasAngles && itAngles != theAngles.end() ) - aPPs.back().SetAngle( *itAngles++ ); - } + for ( size_t i = 0; i < pathNodes.size(); ++i ) + { + ExtrusParam::PathPoint & point = points[ i ]; - TNodeOfNodeListMap mapNewNodes; - TElemOfVecOfNnlmiMap mapElemNewNodes; - TTElemOfElemListMap newElemsMap; - TIDSortedElemSet::iterator itElem; - // source elements for each generated one - SMESH_SequenceOfElemPtr srcElems, srcNodes; + point.myPnt = SMESH_NodeXYZ( pathNodes[ i ]); - // 3. Center of rotation aV0 - gp_Pnt aV0 = theRefPoint; - if ( !theHasRefPoint ) - { - gp_XYZ aGC( 0.,0.,0. ); - TIDSortedElemSet newNodes; + if ( angle != theAngles.end() ) + point.myAngle = *angle++; - for ( int is2ndSet = 0; is2ndSet < 2; ++is2ndSet ) + tangent.SetCoord( 0,0,0 ); + const int shapeID = pathNodes[ i ]->GetShapeID(); + const TopoDS_Shape& shape = pathMeshDS->IndexToShape( shapeID ); + TopAbs_ShapeEnum shapeType = shape.IsNull() ? TopAbs_SHAPE : shape.ShapeType(); + switch ( shapeType ) { - TIDSortedElemSet& theElements = theElemSets[ is2ndSet ]; - itElem = theElements.begin(); - for ( ; itElem != theElements.end(); itElem++ ) + case TopAbs_EDGE: + { + TopoDS_Edge edge = TopoDS::Edge( shape ); + id2factor = edgeID2OriFactor.insert( std::make_pair( shapeID, 0 )).first; + if ( id2factor->second == 0 ) { - const SMDS_MeshElement* elem = *itElem; - SMDS_ElemIteratorPtr itN = elem->nodesIterator(); - while ( itN->more() ) { - const SMDS_MeshElement* node = itN->next(); - if ( newNodes.insert( node ).second ) - aGC += SMESH_NodeXYZ( node ); - } + if ( i ) id2factor->second = getOriFactor( edge, pathNodes[i-1], pathNodes[i], pathHelper ); + else id2factor->second = getOriFactor( edge, pathNodes[i], pathNodes[i+1], pathHelper ); } + double u = pathHelper.GetNodeU( edge, pathNodes[i] ), u0, u1; + Handle(Geom_Curve) curve = BRep_Tool::Curve( edge, u0, u1 ); + curve->D1( u, p, tangent ); + tangent *= id2factor->second; + break; } - aGC /= newNodes.size(); - aV0.SetXYZ( aGC ); - } // if (!theHasRefPoint) { - - // 4. Processing the elements - SMESHDS_Mesh* aMesh = GetMeshDS(); - list emptyList; - - setElemsFirst( theElemSets ); - for ( int is2ndSet = 0; is2ndSet < 2; ++is2ndSet ) - { - TIDSortedElemSet& theElements = theElemSets[ is2ndSet ]; - for ( itElem = theElements.begin(); itElem != theElements.end(); itElem++ ) + case TopAbs_VERTEX: { - const SMDS_MeshElement* elem = *itElem; - - vector & newNodesItVec = mapElemNewNodes[ elem ]; - newNodesItVec.reserve( elem->NbNodes() ); - - // loop on elem nodes - int nodeIndex = -1; - SMDS_ElemIteratorPtr itN = elem->nodesIterator(); - while ( itN->more() ) + int nbEdges = 0; + PShapeIteratorPtr shapeIt = pathHelper.GetAncestors( shape, *theTrackMesh, TopAbs_EDGE ); + while ( const TopoDS_Shape* edgePtr = shapeIt->next() ) { - ++nodeIndex; - // check if a node has been already processed - const SMDS_MeshNode* node = cast2Node( itN->next() ); - TNodeOfNodeListMap::iterator nIt = mapNewNodes.insert( make_pair( node, emptyList )).first; - list& listNewNodes = nIt->second; - if ( listNewNodes.empty() ) + int edgeID = pathMeshDS->ShapeToIndex( *edgePtr ); + for ( int di = -1; di <= 0; ++di ) { - // make new nodes - Standard_Real aAngle1x, aAngleT1T0, aTolAng; - gp_Pnt aP0x, aP1x, aPN0, aPN1, aV0x, aV1x; - gp_Ax1 anAx1, anAxT1T0; - gp_Dir aDT1x, aDT0x, aDT1T0; - - aTolAng=1.e-4; - - aV0x = aV0; - aPN0 = SMESH_NodeXYZ( node ); - - const SMESH_MeshEditor_PathPoint& aPP0 = aPPs[0]; - aP0x = aPP0.Pnt(); - aDT0x= aPP0.Tangent(); - - for ( int j = 1; j < aNbTP; ++j ) { - const SMESH_MeshEditor_PathPoint& aPP1 = aPPs[j]; - aP1x = aPP1.Pnt(); - aDT1x = aPP1.Tangent(); - aAngle1x = aPP1.Angle(); - - gp_Trsf aTrsf, aTrsfRot, aTrsfRotT1T0; - // Translation - gp_Vec aV01x( aP0x, aP1x ); - aTrsf.SetTranslation( aV01x ); - - // translated point - aV1x = aV0x.Transformed( aTrsf ); - aPN1 = aPN0.Transformed( aTrsf ); - - // rotation 1 [ T1,T0 ] - aAngleT1T0=-aDT1x.Angle( aDT0x ); - if (fabs(aAngleT1T0) > aTolAng) + size_t j = i + di; + if ( j < pathEdges.size() && edgeID == pathEdges[ j ]->GetShapeID() ) + { + TopoDS_Edge edge = TopoDS::Edge( *edgePtr ); + id2factor = edgeID2OriFactor.insert( std::make_pair( edgeID, 0 )).first; + if ( id2factor->second == 0 ) { - aDT1T0=aDT1x^aDT0x; - anAxT1T0.SetLocation( aV1x ); - anAxT1T0.SetDirection( aDT1T0 ); - aTrsfRotT1T0.SetRotation( anAxT1T0, aAngleT1T0 ); - - aPN1 = aPN1.Transformed( aTrsfRotT1T0 ); - } - - // rotation 2 - if ( theHasAngles ) { - anAx1.SetLocation( aV1x ); - anAx1.SetDirection( aDT1x ); - aTrsfRot.SetRotation( anAx1, aAngle1x ); - - aPN1 = aPN1.Transformed( aTrsfRot ); + if ( j < i ) + id2factor->second = getOriFactor( edge, pathNodes[i-1], pathNodes[i], pathHelper ); + else + id2factor->second = getOriFactor( edge, pathNodes[i], pathNodes[i+1], pathHelper ); } - - // make new node - if ( elem->IsQuadratic() && !elem->IsMediumNode(node) ) + double u = pathHelper.GetNodeU( edge, pathNodes[i] ), u0, u1; + Handle(Geom_Curve) curve = BRep_Tool::Curve( edge, u0, u1 ); + gp_Vec du; + curve->D1( u, p, du ); + double size2 = du.SquareMagnitude(); + if ( du.SquareMagnitude() > tol2 ) { - // create additional node - gp_XYZ midP = 0.5 * ( aPN1.XYZ() + aPN0.XYZ() ); - const SMDS_MeshNode* newNode = aMesh->AddNode( midP.X(), midP.Y(), midP.Z() ); - myLastCreatedNodes.push_back(newNode); - srcNodes.push_back( node ); - listNewNodes.push_back( newNode ); + tangent += du.Divided( Sqrt( size2 )) * id2factor->second; + nbEdges++; } - const SMDS_MeshNode* newNode = aMesh->AddNode( aPN1.X(), aPN1.Y(), aPN1.Z() ); - myLastCreatedNodes.push_back(newNode); - srcNodes.push_back( node ); - listNewNodes.push_back( newNode ); - - aPN0 = aPN1; - aP0x = aP1x; - aV0x = aV1x; - aDT0x = aDT1x; + break; } } - else if( elem->IsQuadratic() && !elem->IsMediumNode(node) ) + } + if ( nbEdges > 0 ) + break; + } + // fall through + default: + { + for ( int di = -1; di <= 1; di += 2 ) + { + size_t j = i + di; + if ( j < pathNodes.size() ) { - // if current elem is quadratic and current node is not medium - // we have to check - may be it is needed to insert additional nodes - list< const SMDS_MeshNode* > & listNewNodes = nIt->second; - if ((int) listNewNodes.size() == aNbTP-1 ) - { - vector aNodes(2*(aNbTP-1)); - gp_XYZ P(node->X(), node->Y(), node->Z()); - list< const SMDS_MeshNode* >::iterator it = listNewNodes.begin(); - int i; - for(i=0; iX() + P.X() )/2.; - double y = ( N->Y() + P.Y() )/2.; - double z = ( N->Z() + P.Z() )/2.; - const SMDS_MeshNode* newN = aMesh->AddNode(x,y,z); - srcNodes.push_back( node ); - myLastCreatedNodes.push_back(newN); - aNodes[2*i] = newN; - aNodes[2*i+1] = N; - P = gp_XYZ(N->X(),N->Y(),N->Z()); - } - listNewNodes.clear(); - for(i=0; i<2*(aNbTP-1); i++) { - listNewNodes.push_back(aNodes[i]); - } - } + gp_Vec dir( point.myPnt, SMESH_NodeXYZ( pathNodes[ j ])); + double size2 = dir.SquareMagnitude(); + if ( size2 > tol2 ) + tangent += dir.Divided( Sqrt( size2 )) * di; } + } + } + } // switch ( shapeType ) + + if ( tangent.SquareMagnitude() < tol2 ) + return EXTR_CANT_GET_TANGENT; - newNodesItVec.push_back( nIt ); - } + point.myTgt = tangent; - // make new elements - sweepElement( elem, newNodesItVec, newElemsMap[elem], aNbTP-1, srcElems ); - } - } + } // loop on pathNodes - makeWalls( mapNewNodes, newElemsMap, mapElemNewNodes, theElemSets[0], aNbTP-1, srcElems ); - if ( theMakeGroups ) - generateGroups( srcNodes, srcElems, "extruded"); + ExtrusParam nodeMaker( points, theRefPoint, theScales, theMakeGroups ); + TTElemOfElemListMap newElemsMap; + + ExtrusionSweep( theElements, nodeMaker, newElemsMap ); return EXTR_OK; } - //======================================================================= //function : linearAngleVariation //purpose : spread values over nbSteps @@ -6491,6 +6608,33 @@ void SMESH_MeshEditor::linearAngleVariation(const int nbSteps, } } +//======================================================================= +//function : linearScaleVariation +//purpose : spread values over nbSteps +//======================================================================= + +void SMESH_MeshEditor::linearScaleVariation(const int theNbSteps, + std::list& theScales) +{ + int nbScales = theScales.size(); + std::vector myScales; + myScales.reserve( theNbSteps ); + std::list::const_iterator scale = theScales.begin(); + double prevScale = 1.0; + for ( int iSc = 1; scale != theScales.end(); ++scale, ++iSc ) + { + int iStep = int( iSc / double( nbScales ) * theNbSteps + 0.5 ); + int stDelta = Max( 1, iStep - myScales.size()); + double scDelta = ( *scale - prevScale ) / stDelta; + for ( int iStep = 0; iStep < stDelta; ++iStep ) + { + myScales.push_back( prevScale + scDelta ); + prevScale = myScales.back(); + } + prevScale = *scale; + } + theScales.assign( myScales.begin(), myScales.end() ); +} //================================================================================ /*! @@ -6733,7 +6877,7 @@ SMESH_MeshEditor::Transform (TIDSortedElemSet & theElems, * \param [in] theValue - offset value * \param [out] theTgtMesh - a mesh to add offset elements to * \param [in] theMakeGroups - to generate groups - * \return PGroupIDs - IDs of created groups + * \return PGroupIDs - IDs of created groups. NULL means failure */ //================================================================================ @@ -6752,8 +6896,8 @@ SMESH_MeshEditor::PGroupIDs SMESH_MeshEditor::Offset( TIDSortedElemSet & theElem if ( theElements.empty() ) eIt = meshDS->elementsIterator( SMDSAbs_Face ); else eIt = SMESHUtils::elemSetIterator( theElements ); - SMESH_MeshAlgos::TEPairVec new2OldFaces; - SMESH_MeshAlgos::TNPairVec new2OldNodes; + SMESH_MeshAlgos::TElemIntPairVec new2OldFaces; + SMESH_MeshAlgos::TNodeIntPairVec new2OldNodes; std::unique_ptr< SMDS_Mesh > offsetMesh ( SMESH_MeshAlgos::MakeOffset( eIt, *meshDS, theValue, theFixSelfIntersection, @@ -6785,16 +6929,17 @@ SMESH_MeshEditor::PGroupIDs SMESH_MeshEditor::Offset( TIDSortedElemSet & theElem // copy offsetMesh to theTgtMesh - int idShift = meshDS->MaxNodeID(); + smIdType idShift = meshDS->MaxNodeID(); for ( size_t i = 0; i < new2OldNodes.size(); ++i ) if ( const SMDS_MeshNode* n = new2OldNodes[ i ].first ) { - if ( n->NbInverseElements() > 0 ) + + if (!SALOME::VerbosityActivated() || n->NbInverseElements() > 0 ) { const SMDS_MeshNode* n2 = tgtMeshDS->AddNodeWithID( n->X(), n->Y(), n->Z(), idShift + n->GetID() ); myLastCreatedNodes.push_back( n2 ); - srcNodes.push_back( new2OldNodes[ i ].second ); + srcNodes.push_back( meshDS->FindNode( new2OldNodes[ i ].second )); } } @@ -6810,7 +6955,7 @@ SMESH_MeshEditor::PGroupIDs SMESH_MeshEditor::Offset( TIDSortedElemSet & theElem elemType.myNodes.push_back( tgtMeshDS->FindNode( idShift + n2->GetID() )); } tgtEditor.AddElement( elemType.myNodes, elemType ); - srcElems.push_back( new2OldFaces[ i ].second ); + srcElems.push_back( meshDS->FindElement( new2OldFaces[ i ].second )); } myLastCreatedElems.swap( tgtEditor.myLastCreatedElems ); @@ -6818,6 +6963,8 @@ SMESH_MeshEditor::PGroupIDs SMESH_MeshEditor::Offset( TIDSortedElemSet & theElem PGroupIDs newGroupIDs; if ( theMakeGroups ) newGroupIDs = generateGroups( srcNodes, srcElems, "offset", theTgtMesh, false ); + else + newGroupIDs.reset( new std::list< int > ); return newGroupIDs; } @@ -7151,7 +7298,7 @@ void SMESH_MeshEditor::MergeNodes (TListOfListOfNodes & theGroupsOfNodes, TNodeNodeMap nodeNodeMap; // node to replace - new node set elems; // all elements with changed nodes - list< int > rmElemIds, rmNodeIds; + list< smIdType > rmElemIds, rmNodeIds; vector< ElemFeatures > newElemDefs; // Fill nodeNodeMap and elems @@ -7250,6 +7397,7 @@ void SMESH_MeshEditor::MergeNodes (TListOfListOfNodes & theGroupsOfNodes, { const SMDS_MeshElement* elem = *eIt; SMESHDS_SubMesh* sm = mesh->MeshElements( elem->getshapeId() ); + bool marked = elem->isMarked(); bool keepElem = applyMerge( elem, newElemDefs, nodeNodeMap, /*noHoles=*/false ); if ( !keepElem ) @@ -7257,9 +7405,19 @@ void SMESH_MeshEditor::MergeNodes (TListOfListOfNodes & theGroupsOfNodes, for ( size_t i = 0; i < newElemDefs.size(); ++i ) { - if ( i > 0 || !mesh->ChangeElementNodes( elem, - & newElemDefs[i].myNodes[0], - newElemDefs[i].myNodes.size() )) + bool elemChanged = false; + if ( i == 0 ) + { + if ( elem->GetGeomType() == SMDSGeom_POLYHEDRA ) + elemChanged = mesh->ChangePolyhedronNodes( elem, + newElemDefs[i].myNodes, + newElemDefs[i].myPolyhedQuantities ); + else + elemChanged = mesh->ChangeElementNodes( elem, + & newElemDefs[i].myNodes[0], + newElemDefs[i].myNodes.size() ); + } + if ( i > 0 || !elemChanged ) { if ( i == 0 ) { @@ -7276,6 +7434,8 @@ void SMESH_MeshEditor::MergeNodes (TListOfListOfNodes & theGroupsOfNodes, sm->AddElement( newElem ); if ( elem != newElem ) ReplaceElemInGroups( elem, newElem, mesh ); + if ( marked && newElem ) + newElem->setIsMarked( true ); } } } @@ -7426,6 +7586,7 @@ bool SMESH_MeshEditor::applyMerge( const SMDS_MeshElement* elem, // each face has to be analyzed in order to check volume validity if ( const SMDS_MeshVolume* aPolyedre = SMDS_Mesh::DownCast< SMDS_MeshVolume >( elem )) { + toRemove = false; int nbFaces = aPolyedre->NbFaces(); vector& poly_nodes = newElemDefs[0].myNodes; @@ -7746,14 +7907,18 @@ bool SMESH_MeshEditor::applyMerge( const SMDS_MeshElement* elem, // purpose : allow comparing elements basing on their nodes // ======================================================== -class ComparableElement : public boost::container::flat_set< int > +struct ComparableElementHasher; + +class ComparableElement : public boost::container::flat_set< smIdType > { - typedef boost::container::flat_set< int > int_set; + typedef boost::container::flat_set< smIdType > int_set; const SMDS_MeshElement* myElem; - int mySumID; + smIdType mySumID; mutable int myGroupID; + friend ComparableElementHasher; + public: ComparableElement( const SMDS_MeshElement* theElem ): @@ -7762,7 +7927,7 @@ public: this->reserve( theElem->NbNodes() ); for ( SMDS_ElemIteratorPtr nodeIt = theElem->nodesIterator(); nodeIt->more(); ) { - int id = nodeIt->next()->GetID(); + smIdType id = nodeIt->next()->GetID(); mySumID += id; this->insert( id ); } @@ -7774,23 +7939,38 @@ public: //int& GroupID() const { return const_cast< int& >( myGroupID ); } ComparableElement( const ComparableElement& theSource ) // move copy + : int_set() { ComparableElement& src = const_cast< ComparableElement& >( theSource ); - (int_set&) (*this ) = boost::move( src ); + (int_set&) (*this ) = std::move( src ); myElem = src.myElem; mySumID = src.mySumID; myGroupID = src.myGroupID; } +}; +struct ComparableElementHasher +{ +#if OCC_VERSION_LARGE < 0x07080000 static int HashCode(const ComparableElement& se, int limit ) { - return ::HashCode( se.mySumID, limit ); + return ::HashCode( FromSmIdType(se.mySumID), limit ); } static Standard_Boolean IsEqual(const ComparableElement& se1, const ComparableElement& se2 ) { return ( se1 == se2 ); } +#else + size_t operator()(const ComparableElement& se) const + { + return static_cast(FromSmIdType(se.mySumID)); + } + bool operator()(const ComparableElement& se1, const ComparableElement& se2) const + { + return ( se1 == se2 ); + } +#endif }; //======================================================================= @@ -7808,8 +7988,8 @@ void SMESH_MeshEditor::FindEqualElements( TIDSortedElemSet & theElements, if ( theElements.empty() ) elemIt = GetMeshDS()->elementsIterator(); else elemIt = SMESHUtils::elemSetIterator( theElements ); - typedef NCollection_Map< ComparableElement, ComparableElement > TMapOfElements; - typedef std::list TGroupOfElems; + typedef NCollection_Map< ComparableElement, ComparableElementHasher > TMapOfElements; + typedef std::list TGroupOfElems; TMapOfElements mapOfElements; std::vector< TGroupOfElems > arrayOfGroups; TGroupOfElems groupOfElems; @@ -7817,6 +7997,8 @@ void SMESH_MeshEditor::FindEqualElements( TIDSortedElemSet & theElements, while ( elemIt->more() ) { const SMDS_MeshElement* curElem = elemIt->next(); + if ( curElem->IsNull() ) + continue; ComparableElement compElem = curElem; // check uniqueness const ComparableElement& elemInSet = mapOfElements.Added( compElem ); @@ -7853,7 +8035,7 @@ void SMESH_MeshEditor::MergeElements(TListOfListOfElementsID & theGroupsOfElemen { ClearLastCreated(); - typedef list TListOfIDs; + typedef list TListOfIDs; TListOfIDs rmElemIds; // IDs of elems to remove SMESHDS_Mesh* aMesh = GetMeshDS(); @@ -7952,35 +8134,37 @@ bool SMESH_MeshEditor::FindFreeBorder (const SMDS_MeshNode* theFirst theNodes.push_back( theSecondNode ); const SMDS_MeshNode *nIgnore = theFirstNode, *nStart = theSecondNode; - TIDSortedElemSet foundElems; + //TIDSortedElemSet foundElems; bool needTheLast = ( theLastNode != 0 ); + vector nodes; + while ( nStart != theLastNode ) { if ( nStart == theFirstNode ) return !needTheLast; - // find all free border faces sharing form nStart + // find all free border faces sharing nStart list< const SMDS_MeshElement* > curElemList; list< const SMDS_MeshNode* > nStartList; SMDS_ElemIteratorPtr invElemIt = nStart->GetInverseElementIterator(SMDSAbs_Face); while ( invElemIt->more() ) { const SMDS_MeshElement* e = invElemIt->next(); - if ( e == curElem || foundElems.insert( e ).second ) { + //if ( e == curElem || foundElems.insert( e ).second ) // e can encounter twice in border + { // get nodes - int iNode = 0, nbNodes = e->NbNodes(); - vector nodes( nbNodes+1 ); nodes.assign( SMDS_MeshElement::iterator( e->interlacedNodesIterator() ), SMDS_MeshElement::iterator() ); nodes.push_back( nodes[ 0 ]); // check 2 links + int iNode = 0, nbNodes = nodes.size() - 1; for ( iNode = 0; iNode < nbNodes; iNode++ ) - if (((nodes[ iNode ] == nStart && nodes[ iNode + 1] != nIgnore ) || - (nodes[ iNode + 1] == nStart && nodes[ iNode ] != nIgnore )) && - ControlFreeBorder( &nodes[ iNode ], e->GetID() )) + if ((( nodes[ iNode ] == nStart && nodes[ iNode + 1] != nIgnore ) || + ( nodes[ iNode + 1] == nStart && nodes[ iNode ] != nIgnore )) && + ( ControlFreeBorder( &nodes[ iNode ], e->GetID() ))) { - nStartList.push_back( nodes[ iNode + ( nodes[ iNode ] == nStart ? 1 : 0 )]); + nStartList.push_back( nodes[ iNode + ( nodes[ iNode ] == nStart )]); curElemList.push_back( e ); } } @@ -8032,7 +8216,7 @@ bool SMESH_MeshEditor::FindFreeBorder (const SMDS_MeshNode* theFirst else if ( !contNodes[0].empty() && !contNodes[1].empty() ) { // choice: clear a worse one int iLongest = ( contNodes[0].size() < contNodes[1].size() ? 1 : 0 ); - int iWorse = ( needTheLast ? 1 - iLongest : iLongest ); + int iWorse = ( needTheLast ? 1 - iLongest : iLongest ); contNodes[ iWorse ].clear(); contFaces[ iWorse ].clear(); } @@ -8043,13 +8227,11 @@ bool SMESH_MeshEditor::FindFreeBorder (const SMDS_MeshNode* theFirst // push_back the best free border cNL = & contNodes[ contNodes[0].empty() ? 1 : 0 ]; cFL = & contFaces[ contFaces[0].empty() ? 1 : 0 ]; - theNodes.pop_back(); // remove nIgnore + //theNodes.pop_back(); // remove nIgnore theNodes.pop_back(); // remove nStart - theFaces.pop_back(); // remove curElem - list< const SMDS_MeshNode* >::iterator nIt = cNL->begin(); - list< const SMDS_MeshElement* >::iterator fIt = cFL->begin(); - for ( ; nIt != cNL->end(); nIt++ ) theNodes.push_back( *nIt ); - for ( ; fIt != cFL->end(); fIt++ ) theFaces.push_back( *fIt ); + //theFaces.pop_back(); // remove curElem + theNodes.splice( theNodes.end(), *cNL ); + theFaces.splice( theFaces.end(), *cFL ); return true; } // several continuations found @@ -8390,6 +8572,10 @@ SMESH_MeshEditor::SewFreeBorder (const SMDS_MeshNode* theBordFirstNode, nIt[0] = nSide[0].begin(); eIt[0] = eSide[0].begin(); nIt[1] = nSide[1].begin(); eIt[1] = eSide[1].begin(); + // element can be split while iterating on border if it has two edges in the border + std::map< const SMDS_MeshElement* , const SMDS_MeshElement* > elemReplaceMap; + std::map< const SMDS_MeshElement* , const SMDS_MeshElement* >::iterator elemReplaceMapIt; + TElemOfNodeListMap insertMap; TElemOfNodeListMap::iterator insertMapIt; // insertMap is @@ -8437,12 +8623,15 @@ SMESH_MeshEditor::SewFreeBorder (const SMDS_MeshNode* theBordFirstNode, const SMDS_MeshNode* nIns = *nIt [ 1 - intoBord ]; if ( intoBord == 1 ) { // move node of the border to be on a link of elem of the side - gp_XYZ p1 (n1->X(), n1->Y(), n1->Z()); - gp_XYZ p2 (n2->X(), n2->Y(), n2->Z()); + SMESH_NodeXYZ p1( n1 ), p2( n2 ); double ratio = du / ( param[ 1 ][ i[1] ] - param[ 1 ][ i[1]-1 ]); gp_XYZ p = p2 * ( 1 - ratio ) + p1 * ratio; GetMeshDS()->MoveNode( nIns, p.X(), p.Y(), p.Z() ); } + elemReplaceMapIt = elemReplaceMap.find( elem ); + if ( elemReplaceMapIt != elemReplaceMap.end() ) + elem = elemReplaceMapIt->second; + insertMapIt = insertMap.find( elem ); bool notFound = ( insertMapIt == insertMap.end() ); bool otherLink = ( !notFound && (*insertMapIt).second.front() != n1 ); @@ -8465,8 +8654,10 @@ SMESH_MeshEditor::SewFreeBorder (const SMDS_MeshNode* theBordFirstNode, UpdateVolumes(n12, n22, nodeList); } // 3. find an element appeared on n1 and n2 after the insertion - insertMap.erase( elem ); - elem = findAdjacentFace( n1, n2, 0 ); + insertMap.erase( insertMapIt ); + const SMDS_MeshElement* elem2 = findAdjacentFace( n1, n2, 0 ); + elemReplaceMap.insert( std::make_pair( elem, elem2 )); + elem = elem2; } if ( notFound || otherLink ) { // add element and nodes of the side into the insertMap @@ -8499,6 +8690,7 @@ SMESH_MeshEditor::SewFreeBorder (const SMDS_MeshNode* theBordFirstNode, { const SMDS_MeshElement* elem = (*insertMapIt).first; list & nodeList = (*insertMapIt).second; + if ( nodeList.size() < 3 ) continue; const SMDS_MeshNode* n1 = nodeList.front(); nodeList.pop_front(); const SMDS_MeshNode* n2 = nodeList.front(); nodeList.pop_front(); @@ -8558,8 +8750,8 @@ SMESH_MeshEditor::SewFreeBorder (const SMDS_MeshNode* theBordFirstNode, TListOfListOfElementsID::iterator itGroups = equalGroups.begin(); for ( ; itGroups != equalGroups.end(); ++itGroups ) { - list< int >& group = *itGroups; - list< int >::iterator id = group.begin(); + list< smIdType >& group = *itGroups; + list< smIdType >::iterator id = group.begin(); for ( ++id; id != group.end(); ++id ) if ( const SMDS_MeshElement* seg = GetMeshDS()->FindElement( *id )) segments.erase( seg ); @@ -8991,12 +9183,12 @@ namespace */ //======================================================================= -int SMESH_MeshEditor::convertElemToQuadratic(SMESHDS_SubMesh * theSm, - SMESH_MesherHelper& theHelper, - const bool theForce3d) +smIdType SMESH_MeshEditor::convertElemToQuadratic(SMESHDS_SubMesh * theSm, + SMESH_MesherHelper& theHelper, + const bool theForce3d) { //MESSAGE("convertElemToQuadratic"); - int nbElem = 0; + smIdType nbElem = 0; if( !theSm ) return nbElem; vector nbNodeInFaces; @@ -9049,7 +9241,7 @@ int SMESH_MeshEditor::convertElemToQuadratic(SMESHDS_SubMesh * theSm, } // get elem data needed to re-create it // - const int id = elem->GetID(); + const smIdType id = elem->GetID(); const int nbNodes = elem->NbCornerNodes(); nodes.assign(elem->begin_nodes(), elem->end_nodes()); if ( aGeomType == SMDSEntity_Polyhedra ) @@ -9144,7 +9336,7 @@ void SMESH_MeshEditor::ConvertToQuadratic(const bool theForce3d, const bool theT aHelper.ToFixNodeParameters( true ); // convert elements assigned to sub-meshes - int nbCheckedElems = 0; + smIdType nbCheckedElems = 0; if ( myMesh->HasShapeToMesh() ) { if ( SMESH_subMesh *aSubMesh = myMesh->GetSubMeshContaining(myMesh->GetShapeToMesh())) @@ -9161,7 +9353,7 @@ void SMESH_MeshEditor::ConvertToQuadratic(const bool theForce3d, const bool theT } // convert elements NOT assigned to sub-meshes - int totalNbElems = meshDS->NbEdges() + meshDS->NbFaces() + meshDS->NbVolumes(); + smIdType totalNbElems = meshDS->NbEdges() + meshDS->NbFaces() + meshDS->NbVolumes(); if ( nbCheckedElems < totalNbElems ) // not all elements are in sub-meshes { aHelper.SetElementsOnShape(false); @@ -9174,7 +9366,7 @@ void SMESH_MeshEditor::ConvertToQuadratic(const bool theForce3d, const bool theT const SMDS_MeshEdge* edge = aEdgeItr->next(); if ( !edge->IsQuadratic() ) { - int id = edge->GetID(); + smIdType id = edge->GetID(); const SMDS_MeshNode* n1 = edge->GetNode(0); const SMDS_MeshNode* n2 = edge->GetNode(1); @@ -9215,7 +9407,7 @@ void SMESH_MeshEditor::ConvertToQuadratic(const bool theForce3d, const bool theT if ( alreadyOK ) continue; - const int id = face->GetID(); + const smIdType id = face->GetID(); vector nodes ( face->begin_nodes(), face->end_nodes()); meshDS->RemoveFreeElement(face, smDS, /*fromGroups=*/false); @@ -9271,7 +9463,7 @@ void SMESH_MeshEditor::ConvertToQuadratic(const bool theForce3d, const bool theT continue; } } - const int id = volume->GetID(); + const smIdType id = volume->GetID(); vector nodes (volume->begin_nodes(), volume->end_nodes()); if ( type == SMDSEntity_Polyhedra ) nbNodeInFaces = static_cast(volume)->GetQuantities(); @@ -9451,7 +9643,7 @@ void SMESH_MeshEditor::ConvertToQuadratic(const bool theForce3d, if ( alreadyOK ) continue; const SMDSAbs_ElementType type = elem->GetType(); - const int id = elem->GetID(); + const smIdType id = elem->GetID(); const int nbNodes = elem->NbCornerNodes(); vector nodes ( elem->begin_nodes(), elem->end_nodes()); @@ -9512,15 +9704,15 @@ void SMESH_MeshEditor::ConvertToQuadratic(const bool theForce3d, //======================================================================= /*! * \brief Convert quadratic elements to linear ones and remove quadratic nodes - * \return int - nb of checked elements + * \return smIdType - nb of checked elements */ //======================================================================= -int SMESH_MeshEditor::removeQuadElem(SMESHDS_SubMesh * theSm, - SMDS_ElemIteratorPtr theItr, - const int theShapeID) +smIdType SMESH_MeshEditor::removeQuadElem(SMESHDS_SubMesh * theSm, + SMDS_ElemIteratorPtr theItr, + const int /*theShapeID*/) { - int nbElem = 0; + smIdType nbElem = 0; SMESHDS_Mesh* meshDS = GetMeshDS(); ElemFeatures elemType; vector nodes; @@ -9565,7 +9757,7 @@ int SMESH_MeshEditor::removeQuadElem(SMESHDS_SubMesh * theSm, bool SMESH_MeshEditor::ConvertFromQuadratic() { - int nbCheckedElems = 0; + smIdType nbCheckedElems = 0; if ( myMesh->HasShapeToMesh() ) { if ( SMESH_subMesh *aSubMesh = myMesh->GetSubMeshContaining(myMesh->GetShapeToMesh())) @@ -9579,7 +9771,7 @@ bool SMESH_MeshEditor::ConvertFromQuadratic() } } - int totalNbElems = + smIdType totalNbElems = GetMeshDS()->NbEdges() + GetMeshDS()->NbFaces() + GetMeshDS()->NbVolumes(); if ( nbCheckedElems < totalNbElems ) // not all elements are in submeshes { @@ -9618,7 +9810,7 @@ void SMESH_MeshEditor::ConvertFromQuadratic(TIDSortedElemSet& theElements) if ( theElements.empty() ) return; // collect IDs of medium nodes of theElements; some of these nodes will be removed - set mediumNodeIDs; + set mediumNodeIDs; TIDSortedElemSet::iterator eIt = theElements.begin(); for ( ; eIt != theElements.end(); ++eIt ) { @@ -9637,7 +9829,7 @@ void SMESH_MeshEditor::ConvertFromQuadratic(TIDSortedElemSet& theElements) // get remaining medium nodes TIDSortedNodeSet mediumNodes; - set::iterator nIdsIt = mediumNodeIDs.begin(); + set::iterator nIdsIt = mediumNodeIDs.begin(); for ( ; nIdsIt != mediumNodeIDs.end(); ++nIdsIt ) if ( const SMDS_MeshNode* n = GetMeshDS()->FindNode( *nIdsIt )) mediumNodes.insert( mediumNodes.end(), n ); @@ -9665,7 +9857,7 @@ void SMESH_MeshEditor::ConvertFromQuadratic(TIDSortedElemSet& theElements) const SMDS_MeshElement* eComplex = invIt2->next(); if ( eComplex->IsQuadratic() && !allMediumNodesIn( eComplex, mediumNodes)) { - int nbCommonNodes = SMESH_MeshAlgos::GetCommonNodes( e, eComplex ).size(); + int nbCommonNodes = SMESH_MeshAlgos::NbCommonNodes( e, eComplex ); if ( nbCommonNodes == e->NbNodes()) { complexFound = true; @@ -10126,7 +10318,7 @@ SMESH_MeshEditor::SewSideElements (TIDSortedElemSet& theSide1, if ( aResult != SEW_OK) return aResult; - list< int > nodeIDsToRemove; + list< smIdType > nodeIDsToRemove; vector< const SMDS_MeshNode*> nodes; ElemFeatures elemType; @@ -10207,10 +10399,10 @@ SMESH_MeshEditor::FindMatchingNodes(set& theSide1, set * faceSetPtr[] = { &theSide1, &theSide2 }; nReplaceMap.clear(); - if ( theFirstNode1 != theFirstNode2 ) - nReplaceMap.insert( make_pair( theFirstNode1, theFirstNode2 )); - if ( theSecondNode1 != theSecondNode2 ) - nReplaceMap.insert( make_pair( theSecondNode1, theSecondNode2 )); + //if ( theFirstNode1 != theFirstNode2 ) + nReplaceMap.insert( make_pair( theFirstNode1, theFirstNode2 )); + //if ( theSecondNode1 != theSecondNode2 ) + nReplaceMap.insert( make_pair( theSecondNode1, theSecondNode2 )); set< SMESH_TLink > linkSet; // set of nodes where order of nodes is ignored linkSet.insert( SMESH_TLink( theFirstNode1, theSecondNode1 )); @@ -10500,7 +10692,7 @@ namespace // automatically find theAffectedElems for DoubleNodes() if ( SMESH_MeshAlgos::FaceNormal( _elems[1], norm )) avgNorm += norm; - gp_XYZ bordDir( SMESH_NodeXYZ( _nodes[0] ) - SMESH_NodeXYZ( _nodes[1] )); + gp_XYZ bordDir( SMESH_NodeXYZ( this->_nodes[0] ) - SMESH_NodeXYZ( this->_nodes[1] )); norm = bordDir ^ avgNorm; } else @@ -10526,7 +10718,8 @@ namespace // automatically find theAffectedElems for DoubleNodes() if ( maxX < 0 ) { _elems[0]->setIsMarked( false ); - _elems[1]->setIsMarked( true ); + if ( _elems[1] ) + _elems[1]->setIsMarked( true ); } } @@ -10698,7 +10891,11 @@ namespace // automatically find theAffectedElems for DoubleNodes() { fissure.reserve( theElemsOrNodes.size() ); for ( ; elIt != theElemsOrNodes.end(); ++elIt ) + { fissure.push_back( std::move( FissureBorder( *elIt, elemsByFacet ))); + if ( !fissure.back()._elems[1] ) + fissure.pop_back(); + } } if ( fissure.empty() ) return; @@ -10939,7 +11136,7 @@ void SMESH_MeshEditor::DoubleElements( const TIDSortedElemSet& theElements ) } else { - type = (*theElements.begin())->GetType(); + //type = (*theElements.begin())->GetType(); elemIt = SMESHUtils::elemSetIterator( theElements ); } @@ -10954,7 +11151,8 @@ void SMESH_MeshEditor::DoubleElements( const TIDSortedElemSet& theElements ) while ( elemIt->more() ) { const SMDS_MeshElement* elem = elemIt->next(); - if ( elem->GetType() != type || elem->isMarked() ) + if (( type != SMDSAbs_All && elem->GetType() != type ) || + ( elem->isMarked() )) continue; elemType.Init( elem, /*basicOnly=*/false ); @@ -11060,7 +11258,23 @@ bool SMESH_MeshEditor::doubleNodes(SMESHDS_Mesh* theMeshDS, if ( theIsDoubleElem ) AddElement( newNodes, elemType.Init( anElem, /*basicOnly=*/false )); else - theMeshDS->ChangeElementNodes( anElem, &newNodes[ 0 ], newNodes.size() ); + { + // change element nodes + const SMDSAbs_EntityType geomType = anElem->GetEntityType(); + if ( geomType == SMDSEntity_Polyhedra ) + { + // special treatment for polyhedron + const SMDS_MeshVolume* aPolyhedron = SMDS_Mesh::DownCast< SMDS_MeshVolume >( anElem ); + if (!aPolyhedron) { + MESSAGE("Warning: bad volumic element"); + return false; + } + theMeshDS->ChangePolyhedronNodes( anElem, newNodes, aPolyhedron->GetQuantities() ); + } + else + // standard entity type + theMeshDS->ChangeElementNodes( anElem, &newNodes[ 0 ], newNodes.size() ); + } res = true; } @@ -11204,7 +11418,7 @@ namespace { \brief Identify the elements that will be affected by node duplication (actual duplication is not performed). This method is the first step of DoubleNodeElemGroupsInRegion. \param theElems - list of groups of elements (edges or faces) to be replicated - \param theNodesNot - list of groups of nodes not to replicated + \param theNodesNot - list of groups of nodes not to replicate \param theShape - shape to detect affected elements (element which geometric center located on or inside shape). If the shape is null, detection is done on faces orientations (select elements with a gravity center on the side given by faces normals). @@ -11227,8 +11441,8 @@ bool SMESH_MeshEditor::AffectedElemGroupsInRegion( const TIDSortedElemSet& theEl else { const double aTol = Precision::Confusion(); - auto_ptr< BRepClass3d_SolidClassifier> bsc3d; - auto_ptr<_FaceClassifier> aFaceClassifier; + std::unique_ptr< BRepClass3d_SolidClassifier> bsc3d; + std::unique_ptr<_FaceClassifier> aFaceClassifier; if ( theShape.ShapeType() == TopAbs_SOLID ) { bsc3d.reset( new BRepClass3d_SolidClassifier(theShape));; @@ -11346,7 +11560,7 @@ double SMESH_MeshEditor::OrientedAngle(const gp_Pnt& p0, const gp_Pnt& p1, const try { return n2.AngleWithRef(n1, vref); } - catch ( Standard_Failure ) { + catch ( Standard_Failure& ) { } return Max( v1.Magnitude(), v2.Magnitude() ); } @@ -11386,16 +11600,9 @@ bool SMESH_MeshEditor::DoubleNodesOnGroupBoundaries( const std::vector, DownIdCompare> faceDomains; // face --> (id domain --> id volume) - std::mapcelldom; // cell vtkId --> domain + std::map celldom; // cell vtkId --> domain std::map, DownIdCompare> cellDomains; // oldNode --> (id domain --> id cell) std::map > nodeDomains; // oldId --> (domainId --> newId) - faceDomains.clear(); - celldom.clear(); - cellDomains.clear(); - nodeDomains.clear(); - std::map emptyMap; - std::set emptySet; - emptyMap.clear(); //MESSAGE(".. Number of domains :"<GetVtkID(); + vtkIdType vtkId = anElem->GetVtkID(); //MESSAGE(" vtkId " << vtkId << " smdsId " << anElem->GetID()); int neighborsVtkIds[NBMAXNEIGHBORS]; int downIds[NBMAXNEIGHBORS]; @@ -11453,7 +11660,7 @@ bool SMESH_MeshEditor::DoubleNodesOnGroupBoundaries( const std::vectorGetNeighbors(neighborsVtkIds, downIds, downTypes, vtkId); for (int n = 0; n < nbNeighbors; n++) { - int smdsId = meshDS->FromVtkToSmds(neighborsVtkIds[n]); + smIdType smdsId = meshDS->FromVtkToSmds(neighborsVtkIds[n]); const SMDS_MeshElement* elem = meshDS->FindElement(smdsId); if (elem && ! domain.count(elem)) // neighbor is in another domain : face is shared { @@ -11505,14 +11712,13 @@ bool SMESH_MeshEditor::DoubleNodesOnGroupBoundaries( const std::vectorfirst; //MESSAGE(" --- face " << face.cellId); std::set oldNodes; - oldNodes.clear(); grid->GetNodeIds(oldNodes, face.cellId, face.cellType); std::set::iterator itn = oldNodes.begin(); for (; itn != oldNodes.end(); ++itn) { int oldId = *itn; //MESSAGE(" node " << oldId); - vtkCellLinks::Link l = grid->GetCellLinks()->GetLink(oldId); + vtkCellLinks::Link l = (static_cast (grid->GetCellLinks()))->GetLink(oldId); for (int i=0; ifirst; //MESSAGE(" --- face " << face.cellId); std::set oldNodes; - oldNodes.clear(); grid->GetNodeIds(oldNodes, face.cellId, face.cellType); std::set::iterator itn = oldNodes.begin(); for (; itn != oldNodes.end(); ++itn) @@ -11621,7 +11826,6 @@ bool SMESH_MeshEditor::DoubleNodesOnGroupBoundaries( const std::vectorfirst; //MESSAGE(" --- face " << face.cellId); std::set oldNodes; - oldNodes.clear(); grid->GetNodeIds(oldNodes, face.cellId, face.cellType); int nbMultipleNodes = 0; std::set::iterator itn = oldNodes.begin(); @@ -11682,7 +11886,7 @@ bool SMESH_MeshEditor::DoubleNodesOnGroupBoundaries( const std::vectorFromVtkToSmds(vtkVolIds[ivol]); + smIdType smdsId = meshDS->FromVtkToSmds(vtkVolIds[ivol]); const SMDS_MeshElement* elem = meshDS->FindElement(smdsId); if (domain.count(elem)) { @@ -11691,7 +11895,7 @@ bool SMESH_MeshEditor::DoubleNodesOnGroupBoundaries( const std::vectorGetID()); double values[3] = { 0,0,0 }; vtkIdType npts = 0; - vtkIdType* pts = 0; + vtkIdType const *pts(nullptr); grid->GetCellPoints(vtkVolIds[ivol], npts, pts); for ( vtkIdType i = 0; i < npts; ++i ) { @@ -11746,14 +11950,13 @@ bool SMESH_MeshEditor::DoubleNodesOnGroupBoundaries( const std::vector mapOfJunctionGroups; //MESSAGE(".. Creation of elements: simple junction"); - if (createJointElems) + if ( createJointElems ) { - int idg; string joints2DName = "joints2D"; - mapOfJunctionGroups[joints2DName] = this->myMesh->AddGroup(SMDSAbs_Face, joints2DName.c_str(), idg); + mapOfJunctionGroups[joints2DName] = this->myMesh->AddGroup(SMDSAbs_Face, joints2DName.c_str()); SMESHDS_Group *joints2DGrp = dynamic_cast(mapOfJunctionGroups[joints2DName]->GetGroupDS()); string joints3DName = "joints3D"; - mapOfJunctionGroups[joints3DName] = this->myMesh->AddGroup(SMDSAbs_Volume, joints3DName.c_str(), idg); + mapOfJunctionGroups[joints3DName] = this->myMesh->AddGroup(SMDSAbs_Volume, joints3DName.c_str()); SMESHDS_Group *joints3DGrp = dynamic_cast(mapOfJunctionGroups[joints3DName]->GetGroupDS()); itface = faceDomains.begin(); @@ -11762,15 +11965,14 @@ bool SMESH_MeshEditor::DoubleNodesOnGroupBoundaries( const std::vectorfirst; std::set oldNodes; std::set::iterator itn; - oldNodes.clear(); grid->GetNodeIds(oldNodes, face.cellId, face.cellType); - std::map domvol = itface->second; + std::map domvol = itface->second; std::map::iterator itdom = domvol.begin(); - int dom1 = itdom->first; + int dom1 = itdom->first; int vtkVolId = itdom->second; itdom++; - int dom2 = itdom->first; + int dom2 = itdom->first; SMDS_MeshCell *vol = grid->extrudeVolumeFromFace(vtkVolId, dom1, dom2, oldNodes, nodeDomains, nodeQuadDomains); stringstream grpname; @@ -11781,7 +11983,7 @@ bool SMESH_MeshEditor::DoubleNodesOnGroupBoundaries( const std::vectormyMesh->AddGroup(vol->GetType(), namegrp.c_str(), idg); + mapOfJunctionGroups[namegrp] = this->myMesh->AddGroup(vol->GetType(), namegrp.c_str()); SMESHDS_Group *sgrp = dynamic_cast(mapOfJunctionGroups[namegrp]->GetGroupDS()); if (sgrp) sgrp->Add(vol->GetID()); @@ -11814,10 +12016,9 @@ bool SMESH_MeshEditor::DoubleNodesOnGroupBoundaries( const std::vectormyMesh->AddGroup(SMDSAbs_Face, namegrp.c_str(), idg); + mapOfJunctionGroups[namegrp] = this->myMesh->AddGroup(SMDSAbs_Face, namegrp.c_str()); SMESHDS_Group *sgrp = dynamic_cast(mapOfJunctionGroups[namegrp]->GetGroupDS()); if (sgrp) sgrp->Add(face->GetID()); @@ -11828,7 +12029,7 @@ bool SMESH_MeshEditor::DoubleNodesOnGroupBoundaries( const std::vector, std::vector >::iterator ite = edgesMultiDomains.begin(); for (; ite != edgesMultiDomains.end(); ++ite) { - vector nodes = ite->first; + vector nodes = ite->first; vector orderDom = ite->second; vector orderedNodes; if (nodes.size() == 2) @@ -11843,10 +12044,9 @@ bool SMESH_MeshEditor::DoubleNodesOnGroupBoundaries( const std::vectorGetMeshDS()->AddVolumeFromVtkIds(orderedNodes); - int idg; string namegrp = "jointsMultiples"; if (!mapOfJunctionGroups.count(namegrp)) - mapOfJunctionGroups[namegrp] = this->myMesh->AddGroup(SMDSAbs_Volume, namegrp.c_str(), idg); + mapOfJunctionGroups[namegrp] = this->myMesh->AddGroup(SMDSAbs_Volume, namegrp.c_str()); SMESHDS_Group *sgrp = dynamic_cast(mapOfJunctionGroups[namegrp]->GetGroupDS()); if (sgrp) sgrp->Add(vol->GetID()); @@ -11866,10 +12066,9 @@ bool SMESH_MeshEditor::DoubleNodesOnGroupBoundaries( const std::vector, DownIdCompare> faceOrEdgeDom; // cellToModify --> (id domain --> id cell) std::map feDom; // vtk id of cell to modify --> id domain - faceOrEdgeDom.clear(); - feDom.clear(); //MESSAGE(".. Modification of elements"); + SMDSAbs_ElementType domainType = (*theElems[0].begin())->GetType(); for (int idomain = idom0; idomain < nbDomains; idomain++) { std::map >::const_iterator itnod = nodeDomains.begin(); @@ -11877,7 +12076,7 @@ bool SMESH_MeshEditor::DoubleNodesOnGroupBoundaries( const std::vectorfirst; //MESSAGE(" node " << oldId); - vtkCellLinks::Link l = grid->GetCellLinks()->GetLink(oldId); + vtkCellLinks::Link l = (static_cast< vtkCellLinks *>(grid->GetCellLinks()))->GetLink(oldId); for (int i = 0; i < l.ncells; i++) { int vtkId = l.cells[i]; @@ -11887,13 +12086,29 @@ bool SMESH_MeshEditor::DoubleNodesOnGroupBoundaries( const std::vectorGetParentVolumes(volParents, vtkId); + int nbvol = 0; + nbvol = grid->GetParentVolumes(volParents, vtkId); + if ( domainType == SMDSAbs_Volume ) + { + nbvol = grid->GetParentVolumes(volParents, vtkId); + } + else // domainType == SMDSAbs_Face + { + const int nbFaces = grid->getDownArray(vtkType)->getNumberOfUpCells(downId); + const int *upCells = grid->getDownArray(vtkType)->getUpCells(downId); + const unsigned char* upTypes = grid->getDownArray(vtkType)->getUpTypes(downId); + for (int i=0; i< nbFaces; i++) + { + int vtkFaceId = grid->getDownArray( upTypes[i] )->getVtkCellId(upCells[i]); + if (vtkFaceId >= 0) + volParents[nbvol++] = vtkFaceId; + } + } for (int j = 0; j < nbvol; j++) if (celldom.count(volParents[j]) && (celldom[volParents[j]] == idomain)) if (!feDom.count(vtkId)) { feDom[vtkId] = idomain; - faceOrEdgeDom[aCell] = emptyMap; faceOrEdgeDom[aCell][idomain] = vtkId; // affect face or edge to the first domain only //MESSAGE("affect cell " << this->GetMeshDS()->FromVtkToSmds(vtkId) << " domain " << idomain // << " type " << vtkType << " downId " << downId); @@ -11916,7 +12131,6 @@ bool SMESH_MeshEditor::DoubleNodesOnGroupBoundaries( const std::vectorfirst; std::set oldNodes; std::set::iterator itn; - oldNodes.clear(); grid->GetNodeIds(oldNodes, face.cellId, face.cellType); //MESSAGE("examine cell, downId " << face.cellId << " type " << int(face.cellType)); std::map localClonedNodeIds; @@ -11983,10 +12197,7 @@ bool SMESH_MeshEditor::CreateFlatElementsOnFacesGroups(const std::vector clonedNodes; std::map intermediateNodes; - clonedNodes.clear(); - intermediateNodes.clear(); std::map mapOfJunctionGroups; - mapOfJunctionGroups.clear(); for ( size_t idom = 0; idom < theElems.size(); idom++ ) { @@ -12104,10 +12315,9 @@ bool SMESH_MeshEditor::CreateFlatElementsOnFacesGroups(const std::vectormyMesh->AddGroup(SMDSAbs_Volume, namegrp.c_str(), idg); + mapOfJunctionGroups[namegrp] = this->myMesh->AddGroup(SMDSAbs_Volume, namegrp.c_str()); SMESHDS_Group *sgrp = dynamic_cast(mapOfJunctionGroups[namegrp]->GetGroupDS()); if (sgrp) sgrp->Add(vol->GetID()); @@ -12178,10 +12388,10 @@ void SMESH_MeshEditor::CreateHoleSkin(double radius, // --- define groups to build - int idg; // --- group of SMDS volumes + // --- group of SMDS volumes string grpvName = groupName; grpvName += "_vol"; - SMESH_Group *grp = this->myMesh->AddGroup(SMDSAbs_Volume, grpvName.c_str(), idg); + SMESH_Group *grp = this->myMesh->AddGroup(SMDSAbs_Volume, grpvName.c_str()); if (!grp) { MESSAGE("group not created " << grpvName); @@ -12189,10 +12399,10 @@ void SMESH_MeshEditor::CreateHoleSkin(double radius, } SMESHDS_Group *sgrp = dynamic_cast(grp->GetGroupDS()); - int idgs; // --- group of SMDS faces on the skin + // --- group of SMDS faces on the skin string grpsName = groupName; grpsName += "_skin"; - SMESH_Group *grps = this->myMesh->AddGroup(SMDSAbs_Face, grpsName.c_str(), idgs); + SMESH_Group *grps = this->myMesh->AddGroup(SMDSAbs_Face, grpsName.c_str()); if (!grps) { MESSAGE("group not created " << grpsName); @@ -12200,10 +12410,10 @@ void SMESH_MeshEditor::CreateHoleSkin(double radius, } SMESHDS_Group *sgrps = dynamic_cast(grps->GetGroupDS()); - int idgi; // --- group of SMDS faces internal (several shapes) + // --- group of SMDS faces internal (several shapes) string grpiName = groupName; grpiName += "_internalFaces"; - SMESH_Group *grpi = this->myMesh->AddGroup(SMDSAbs_Face, grpiName.c_str(), idgi); + SMESH_Group *grpi = this->myMesh->AddGroup(SMDSAbs_Face, grpiName.c_str()); if (!grpi) { MESSAGE("group not created " << grpiName); @@ -12211,10 +12421,10 @@ void SMESH_MeshEditor::CreateHoleSkin(double radius, } SMESHDS_Group *sgrpi = dynamic_cast(grpi->GetGroupDS()); - int idgei; // --- group of SMDS faces internal (several shapes) + // --- group of SMDS faces internal (several shapes) string grpeiName = groupName; grpeiName += "_internalEdges"; - SMESH_Group *grpei = this->myMesh->AddGroup(SMDSAbs_Edge, grpeiName.c_str(), idgei); + SMESH_Group *grpei = this->myMesh->AddGroup(SMDSAbs_Edge, grpeiName.c_str()); if (!grpei) { MESSAGE("group not created " << grpeiName); @@ -12234,7 +12444,6 @@ void SMESH_MeshEditor::CreateHoleSkin(double radius, std::set setOfVolToCheck; std::vector gpnts; - gpnts.clear(); if (isNodeGroup) // --- a group of nodes is provided : find all the volumes using one or more of this nodes { @@ -12321,7 +12530,6 @@ void SMESH_MeshEditor::CreateHoleSkin(double radius, // Fill the group of inside volumes std::map mapOfNodeDistance2; - mapOfNodeDistance2.clear(); std::set setOfOutsideVol; while (!setOfVolToCheck.empty()) { @@ -12330,7 +12538,7 @@ void SMESH_MeshEditor::CreateHoleSkin(double radius, //MESSAGE("volume to check, vtkId " << vtkId << " smdsId " << meshDS->FromVtkToSmds(vtkId)); bool volInside = false; vtkIdType npts = 0; - vtkIdType* pts = 0; + vtkIdType const *pts(nullptr); grid->GetCellPoints(vtkId, npts, pts); for (int i=0; i > shapeIdToVtkIdSet; // shapeId --> set of vtkId on skin - std::set emptySet; - emptySet.clear(); - std::set shapeIds; + std::set shapeIds; SMDS_ElemIteratorPtr itelem = sgrps->GetElements(); while (itelem->more()) @@ -12510,7 +12716,6 @@ void SMESH_MeshEditor::CreateHoleSkin(double radius, int vtkId = elem->GetVtkID(); if (!shapeIdToVtkIdSet.count(shapeId)) { - shapeIdToVtkIdSet[shapeId] = emptySet; shapeIds.insert(shapeId); } shapeIdToVtkIdSet[shapeId].insert(vtkId); @@ -12518,7 +12723,6 @@ void SMESH_MeshEditor::CreateHoleSkin(double radius, std::map > shapeIdToEdges; // shapeId --> set of downward edges std::set emptyEdges; - emptyEdges.clear(); std::map >::iterator itShape = shapeIdToVtkIdSet.begin(); for (; itShape != shapeIdToVtkIdSet.end(); ++itShape) @@ -12542,7 +12746,7 @@ void SMESH_MeshEditor::CreateHoleSkin(double radius, { if (neighborsVtkIds[n]<0) // only smds faces are considered as neighbors here continue; - int smdsId = meshDS->FromVtkToSmds(neighborsVtkIds[n]); + smIdType smdsId = meshDS->FromVtkToSmds(neighborsVtkIds[n]); const SMDS_MeshElement* elem = meshDS->FindElement(smdsId); if ( shapeIds.count(elem->getshapeId()) && !sgrps->Contains(elem)) // edge : neighbor in the set of shape, not in the group { @@ -12561,7 +12765,6 @@ void SMESH_MeshEditor::CreateHoleSkin(double radius, } std::list order; - order.clear(); if (nodesEdges.size() > 0) { order.push_back(nodesEdges[0]); //MESSAGE(" --- back " << order.back()+1); // SMDS id = VTK id + 1; @@ -12742,7 +12945,8 @@ int SMESH_MeshEditor::MakeBoundaryMesh(const TIDSortedElemSet& elements, bool toCopyElements/*=false*/, bool toCopyExistingBoundary/*=false*/, bool toAddExistingBondary/*= false*/, - bool aroundElements/*= false*/) + bool aroundElements/*= false*/, + bool toCreateAllElements/*= false*/) { SMDSAbs_ElementType missType = (dimension == BND_2DFROM3D) ? SMDSAbs_Face : SMDSAbs_Edge; SMDSAbs_ElementType elemType = (dimension == BND_1DFROM2D) ? SMDSAbs_Face : SMDSAbs_Volume; @@ -12761,7 +12965,6 @@ int SMESH_MeshEditor::MakeBoundaryMesh(const TIDSortedElemSet& elements, SMESH_MeshEditor* presentEditor; SMESH_MeshEditor tgtEditor2( tgtEditor.GetMesh() ); presentEditor = toAddExistingBondary ? &tgtEditor : &tgtEditor2; - SMESH_MesherHelper helper( *myMesh ); const TopAbs_ShapeEnum missShapeType = ( missType==SMDSAbs_Face ? TopAbs_FACE : TopAbs_EDGE ); SMDS_VolumeTool vTool; @@ -12799,8 +13002,9 @@ int SMESH_MeshEditor::MakeBoundaryMesh(const TIDSortedElemSet& elements, const SMDS_MeshElement* otherVol = 0; for ( int iface = 0, n = vTool.NbFaces(); iface < n; iface++ ) { - if ( !vTool.IsFreeFace(iface, &otherVol) && - ( !aroundElements || elements.count( otherVol ))) + if ( !toCreateAllElements && + !vTool.IsFreeFace(iface, &otherVol) && + ( !aroundElements || elements.count( otherVol ))) continue; freeFacets.push_back( iface ); } @@ -12897,8 +13101,8 @@ int SMESH_MeshEditor::MakeBoundaryMesh(const TIDSortedElemSet& elements, for ( inode = 0; inode < srcNodes.size(); ++inode ) tgtNodes[inode] = getNodeWithSameID( tgtMeshDS, srcNodes[inode] ); if ( /*aroundElements && */tgtEditor.GetMeshDS()->FindElement( tgtNodes, - missType, - /*noMedium=*/false)) + missType, + /*noMedium=*/false)) continue; tgtEditor.AddElement( tgtNodes, elemKind.SetPoly( tgtNodes.size()/(iQuad+1) > 4 )); ++nbAddedBnd; @@ -12908,8 +13112,8 @@ int SMESH_MeshEditor::MakeBoundaryMesh(const TIDSortedElemSet& elements, { TConnectivity& nodes = missingBndElems[ i ]; if ( /*aroundElements && */tgtEditor.GetMeshDS()->FindElement( nodes, - missType, - /*noMedium=*/false)) + missType, + /*noMedium=*/false)) continue; SMDS_MeshElement* newElem = tgtEditor.AddElement( nodes, elemKind.SetPoly( nodes.size()/(iQuad+1) > 4 )); @@ -13038,482 +13242,3 @@ void SMESH_MeshEditor::copyPosition( const SMDS_MeshNode* from, default:; } } - -namespace // utils for MakePolyLine -{ - //================================================================================ - /*! - * \brief Sequence of found points and a current point data - */ - struct Path - { - std::vector< gp_XYZ > myPoints; - double myLength; - - int mySrcPntInd; //!< start point index - const SMDS_MeshElement* myFace; - SMESH_NodeXYZ myNode1; - SMESH_NodeXYZ myNode2; - int myNodeInd1; - int myNodeInd2; - double myDot1; - double myDot2; - TIDSortedElemSet myElemSet, myAvoidSet; - - Path(): myLength(0.0), myFace(0) {} - - bool SetCutAtCorner( const SMESH_NodeXYZ& cornerNode, - const SMDS_MeshElement* face, - const gp_XYZ& plnNorm, - const gp_XYZ& plnOrig ); - - void AddPoint( const gp_XYZ& p ); - - bool Extend( const gp_XYZ& plnNorm, const gp_XYZ& plnOrig ); - - bool ReachSamePoint( const Path& other ); - - static void Remove( std::vector< Path > & paths, size_t& i ); - }; - - //================================================================================ - /*! - * \brief Return true if this Path meats another - */ - //================================================================================ - - bool Path::ReachSamePoint( const Path& other ) - { - return ( mySrcPntInd != other.mySrcPntInd && - myFace == other.myFace ); - } - - //================================================================================ - /*! - * \brief Remove a path from a vector - */ - //================================================================================ - - void Path::Remove( std::vector< Path > & paths, size_t& i ) - { - if ( paths.size() > 1 ) - { - size_t j = paths.size() - 1; // last item to be removed - if ( i < j ) - { - paths[ i ].myPoints.swap( paths[ j ].myPoints ); - paths[ i ].myLength = paths[ j ].myLength; - paths[ i ].mySrcPntInd = paths[ j ].mySrcPntInd; - paths[ i ].myFace = paths[ j ].myFace; - paths[ i ].myNode1 = paths[ j ].myNode1; - paths[ i ].myNode2 = paths[ j ].myNode2; - paths[ i ].myNodeInd1 = paths[ j ].myNodeInd1; - paths[ i ].myNodeInd2 = paths[ j ].myNodeInd2; - paths[ i ].myDot1 = paths[ j ].myDot1; - paths[ i ].myDot2 = paths[ j ].myDot2; - } - } - paths.pop_back(); - if ( i > 0 ) - --i; - } - - //================================================================================ - /*! - * \brief Store a point that is at a node of a face if the face is intersected by plane. - * Return false if the node is a sole intersection point of the face and the plane - */ - //================================================================================ - - bool Path::SetCutAtCorner( const SMESH_NodeXYZ& cornerNode, - const SMDS_MeshElement* face, - const gp_XYZ& plnNorm, - const gp_XYZ& plnOrig ) - { - if ( face == myFace ) - return false; - myNodeInd1 = face->GetNodeIndex( cornerNode._node ); - myNodeInd2 = ( myNodeInd1 + 1 ) % face->NbCornerNodes(); - int ind3 = ( myNodeInd1 + 2 ) % face->NbCornerNodes(); - myNode1.Set( face->GetNode( ind3 )); - myNode2.Set( face->GetNode( myNodeInd2 )); - - myDot1 = plnNorm * ( myNode1 - plnOrig ); - myDot2 = plnNorm * ( myNode2 - plnOrig ); - - bool ok = ( myDot1 * myDot2 < 0 ); - if ( !ok && myDot1 * myDot2 == 0 ) - { - ok = ( myDot1 != myDot2 ); - if ( ok && myFace ) - ok = ( myFace->GetNodeIndex(( myDot1 == 0 ? myNode1 : myNode2 )._node ) < 0 ); - } - if ( ok ) - { - myFace = face; - myDot1 = 0; - AddPoint( cornerNode ); - } - return ok; - } - - //================================================================================ - /*! - * \brief Store a point and update myLength - */ - //================================================================================ - - void Path::AddPoint( const gp_XYZ& p ) - { - if ( !myPoints.empty() ) - myLength += ( p - myPoints.back() ).Modulus(); - else - myLength = 0; - myPoints.push_back( p ); - } - - //================================================================================ - /*! - * \brief Try to find the next point - * \param [in] plnNorm - cutting plane normal - * \param [in] plnOrig - cutting plane origin - */ - //================================================================================ - - bool Path::Extend( const gp_XYZ& plnNorm, const gp_XYZ& plnOrig ) - { - int nodeInd3 = ( myNodeInd1 + 1 ) % myFace->NbCornerNodes(); - if ( myNodeInd2 == nodeInd3 ) - nodeInd3 = ( myNodeInd1 + 2 ) % myFace->NbCornerNodes(); - - SMESH_NodeXYZ node3 = myFace->GetNode( nodeInd3 ); - double dot3 = plnNorm * ( node3 - plnOrig ); - - if ( dot3 * myDot1 < 0. ) - { - myNode2 = node3; - myNodeInd2 = nodeInd3; - myDot2 = dot3; - } - else if ( dot3 * myDot2 < 0. ) - { - myNode1 = node3; - myNodeInd1 = nodeInd3; - myDot1 = dot3; - } - else if ( dot3 == 0. ) - { - SMDS_ElemIteratorPtr fIt = node3._node->GetInverseElementIterator(SMDSAbs_Face); - while ( fIt->more() ) - if ( SetCutAtCorner( node3, fIt->next(), plnNorm, plnOrig )) - return true; - return false; - } - else if ( myDot2 == 0. ) - { - SMESH_NodeXYZ node2 = myNode2; // copy as myNode2 changes in SetCutAtCorner() - SMDS_ElemIteratorPtr fIt = node2._node->GetInverseElementIterator(SMDSAbs_Face); - while ( fIt->more() ) - if ( SetCutAtCorner( node2, fIt->next(), plnNorm, plnOrig )) - return true; - return false; - } - - double r = Abs( myDot1 / ( myDot2 - myDot1 )); - AddPoint( myNode1 * ( 1 - r ) + myNode2 * r ); - - myAvoidSet.clear(); - myAvoidSet.insert( myFace ); - myFace = SMESH_MeshAlgos::FindFaceInSet( myNode1._node, myNode2._node, - myElemSet, myAvoidSet, - &myNodeInd1, &myNodeInd2 ); - return myFace; - } - - //================================================================================ - /*! - * \brief Compute a path between two points of PolySegment - */ - struct PolyPathCompute - { - SMESH_MeshEditor::TListOfPolySegments& mySegments; //!< inout PolySegment's - std::vector< Path >& myPaths; //!< path of each of segments to compute - SMESH_Mesh* myMesh; - mutable std::vector< std::string > myErrors; - - PolyPathCompute( SMESH_MeshEditor::TListOfPolySegments& theSegments, - std::vector< Path >& thePaths, - SMESH_Mesh* theMesh): - mySegments( theSegments ), - myPaths( thePaths ), - myMesh( theMesh ), - myErrors( theSegments.size() ) - { - } -#undef SMESH_CAUGHT -#define SMESH_CAUGHT myErrors[i] = - void operator() ( const int i ) const - { - SMESH_TRY; - const_cast< PolyPathCompute* >( this )->Compute( i ); - SMESH_CATCH( SMESH::returnError ); - } -#undef SMESH_CAUGHT - //================================================================================ - /*! - * \brief Compute a path of a given segment - */ - //================================================================================ - - void Compute( const int iSeg ) - { - SMESH_MeshEditor::PolySegment& polySeg = mySegments[ iSeg ]; - - // get a cutting plane - - gp_XYZ p1 = SMESH_NodeXYZ( polySeg.myNode1[0] ); - gp_XYZ p2 = SMESH_NodeXYZ( polySeg.myNode1[1] ); - if ( polySeg.myNode2[0] ) p1 = 0.5 * ( p1 + SMESH_NodeXYZ( polySeg.myNode2[0] )); - if ( polySeg.myNode2[1] ) p2 = 0.5 * ( p2 + SMESH_NodeXYZ( polySeg.myNode2[1] )); - - gp_XYZ plnNorm = ( p1 - p2 ) ^ polySeg.myVector.XYZ(); - gp_XYZ plnOrig = p2; - - // find paths connecting the 2 end points of polySeg - - std::vector< Path > paths; paths.reserve(10); - - // initialize paths - - for ( int iP = 0; iP < 2; ++iP ) // loop on the polySeg end points - { - Path path; - path.mySrcPntInd = iP; - size_t nbPaths = paths.size(); - - if ( polySeg.myNode2[ iP ] && polySeg.myNode2[ iP ] != polySeg.myNode1[ iP ] ) - { - while (( path.myFace = SMESH_MeshAlgos::FindFaceInSet( polySeg.myNode1[ iP ], - polySeg.myNode2[ iP ], - path.myElemSet, - path.myAvoidSet, - &path.myNodeInd1, - &path.myNodeInd2 ))) - { - path.myNode1.Set( polySeg.myNode1[ iP ]); - path.myNode2.Set( polySeg.myNode2[ iP ]); - path.myDot1 = plnNorm * ( path.myNode1 - plnOrig ); - path.myDot2 = plnNorm * ( path.myNode2 - plnOrig ); - path.myPoints.clear(); - path.AddPoint( 0.5 * ( path.myNode1 + path.myNode2 )); - path.myAvoidSet.insert( path.myFace ); - paths.push_back( path ); - } - if ( nbPaths == paths.size() ) - throw SALOME_Exception ( SMESH_Comment("No face edge found by point ") << iP+1 - << " in a PolySegment " << iSeg ); - } - else // an end point is at node - { - std::set nodes; - SMDS_ElemIteratorPtr fIt = polySeg.myNode1[ iP ]->GetInverseElementIterator(SMDSAbs_Face); - while ( fIt->more() ) - { - path.myPoints.clear(); - if ( path.SetCutAtCorner( polySeg.myNode1[ iP ], fIt->next(), plnNorm, plnOrig )) - { - if (( path.myDot1 * path.myDot2 != 0 ) || - ( nodes.insert( path.myDot1 == 0 ? path.myNode1._node : path.myNode2._node ).second )) - paths.push_back( path ); - } - } - } - - // look for a one-segment path - for ( size_t i = 0; i < nbPaths; ++i ) - for ( size_t j = nbPaths; j < paths.size(); ++j ) - if ( paths[i].myFace == paths[j].myFace ) - { - myPaths[ iSeg ].myPoints.push_back( paths[i].myPoints[0] ); - myPaths[ iSeg ].myPoints.push_back( paths[j].myPoints[0] ); - paths.clear(); - } - } - - // extend paths - - myPaths[ iSeg ].myLength = 1e100; - - while ( paths.size() >= 2 ) - { - for ( size_t i = 0; i < paths.size(); ++i ) - { - Path& path = paths[ i ]; - if ( !path.Extend( plnNorm, plnOrig ) || // path reached a mesh boundary - path.myLength > myPaths[ iSeg ].myLength ) // path is longer than others - { - Path::Remove( paths, i ); - continue; - } - - // join paths that reach same point - for ( size_t j = 0; j < paths.size(); ++j ) - { - if ( i != j && paths[i].ReachSamePoint( paths[j] )) - { - double distLast = ( paths[i].myPoints.back() - paths[j].myPoints.back() ).Modulus(); - double fullLength = ( paths[i].myLength + paths[j].myLength + distLast ); - if ( fullLength < myPaths[ iSeg ].myLength ) - { - myPaths[ iSeg ].myLength = fullLength; - std::vector< gp_XYZ > & allPoints = myPaths[ iSeg ].myPoints; - allPoints.swap( paths[i].myPoints ); - allPoints.insert( allPoints.end(), - paths[j].myPoints.rbegin(), - paths[j].myPoints.rend() ); - } - Path::Remove( paths, i ); - Path::Remove( paths, j ); - } - } - } - if ( !paths.empty() && (int) paths[0].myPoints.size() > myMesh->NbFaces() ) - throw SALOME_Exception(LOCALIZED( "Infinite loop in MakePolyLine()")); - } - - if ( myPaths[ iSeg ].myPoints.empty() ) - throw SALOME_Exception( SMESH_Comment("Can't find a full path for PolySegment #") << iSeg ); - - } // PolyPathCompute::Compute() - - }; // struct PolyPathCompute - -} // namespace - -//======================================================================= -//function : MakePolyLine -//purpose : Create a polyline consisting of 1D mesh elements each lying on a 2D element of -// the initial mesh -//======================================================================= - -void SMESH_MeshEditor::MakePolyLine( TListOfPolySegments& theSegments, - SMESHDS_Group* theGroup, - SMESH_ElementSearcher* theSearcher) -{ - std::vector< Path > segPaths( theSegments.size() ); // path of each of segments - - SMESH_ElementSearcher* searcher = theSearcher; - SMESHUtils::Deleter delSearcher; - if ( !searcher ) - { - searcher = SMESH_MeshAlgos::GetElementSearcher( *GetMeshDS() ); - delSearcher._obj = searcher; - } - - // get cutting planes - - std::vector< bool > isVectorOK( theSegments.size(), true ); - const double planarCoef = 0.333; // plane height in planar case - - for ( size_t iSeg = 0; iSeg < theSegments.size(); ++iSeg ) - { - PolySegment& polySeg = theSegments[ iSeg ]; - - gp_XYZ p1 = SMESH_NodeXYZ( polySeg.myNode1[0] ); - gp_XYZ p2 = SMESH_NodeXYZ( polySeg.myNode1[1] ); - if ( polySeg.myNode2[0] ) p1 = 0.5 * ( p1 + SMESH_NodeXYZ( polySeg.myNode2[0] )); - if ( polySeg.myNode2[1] ) p2 = 0.5 * ( p2 + SMESH_NodeXYZ( polySeg.myNode2[1] )); - - gp_XYZ plnNorm = ( p1 - p2 ) ^ polySeg.myVector.XYZ(); - - isVectorOK[ iSeg ] = ( plnNorm.Modulus() > std::numeric_limits::min() ); - if ( !isVectorOK[ iSeg ]) - { - gp_XYZ pMid = 0.5 * ( p1 + p2 ); - const SMDS_MeshElement* face; - polySeg.myMidProjPoint = searcher->Project( pMid, SMDSAbs_Face, &face ); - polySeg.myVector = polySeg.myMidProjPoint.XYZ() - pMid; - - gp_XYZ faceNorm; - SMESH_MeshAlgos::FaceNormal( face, faceNorm ); - - if ( polySeg.myVector.Magnitude() < Precision::Confusion() || - polySeg.myVector * faceNorm < Precision::Confusion() ) - { - polySeg.myVector = faceNorm; - polySeg.myMidProjPoint = pMid + faceNorm * ( p1 - p2 ).Modulus() * planarCoef; - } - } - else - { - polySeg.myVector = plnNorm ^ ( p1 - p2 ); - } - } - - // assure that inverse elements are constructed, avoid their concurrent building in threads - GetMeshDS()->nodesIterator()->next()->NbInverseElements(); - - // find paths - - PolyPathCompute algo( theSegments, segPaths, myMesh ); - OSD_Parallel::For( 0, theSegments.size(), algo, theSegments.size() == 1 ); - - for ( size_t iSeg = 0; iSeg < theSegments.size(); ++iSeg ) - if ( !algo.myErrors[ iSeg ].empty() ) - throw SALOME_Exception( algo.myErrors[ iSeg ].c_str() ); - - // create an 1D mesh - - const SMDS_MeshNode *n, *nPrev = 0; - SMESHDS_Mesh* mesh = GetMeshDS(); - - for ( size_t iSeg = 0; iSeg < theSegments.size(); ++iSeg ) - { - const Path& path = segPaths[iSeg]; - if ( path.myPoints.size() < 2 ) - continue; - - double tol = path.myLength / path.myPoints.size() / 1000.; - if ( !nPrev || ( SMESH_NodeXYZ( nPrev ) - path.myPoints[0] ).SquareModulus() > tol*tol ) - { - nPrev = mesh->AddNode( path.myPoints[0].X(), path.myPoints[0].Y(), path.myPoints[0].Z() ); - myLastCreatedNodes.push_back( nPrev ); - } - for ( size_t iP = 1; iP < path.myPoints.size(); ++iP ) - { - n = mesh->AddNode( path.myPoints[iP].X(), path.myPoints[iP].Y(), path.myPoints[iP].Z() ); - myLastCreatedNodes.push_back( n ); - - const SMDS_MeshElement* elem = mesh->AddEdge( nPrev, n ); - myLastCreatedElems.push_back( elem ); - if ( theGroup ) - theGroup->Add( elem ); - - nPrev = n; - } - - // return a vector - - gp_XYZ pMid = 0.5 * ( path.myPoints[0] + path.myPoints.back() ); - if ( isVectorOK[ iSeg ]) - { - // find the most distance point of a path - double maxDist = 0; - for ( size_t iP = 1; iP < path.myPoints.size(); ++iP ) - { - double dist = Abs( theSegments[iSeg].myVector * ( path.myPoints[iP] - path.myPoints[0] )); - if ( dist > maxDist ) - { - maxDist = dist; - theSegments[iSeg].myMidProjPoint = path.myPoints[iP]; - } - } - if ( maxDist < Precision::Confusion() ) // planar case - theSegments[iSeg].myMidProjPoint = - pMid + theSegments[iSeg].myVector.XYZ().Normalized() * path.myLength * planarCoef; - } - theSegments[iSeg].myVector = gp_Vec( pMid, theSegments[iSeg].myMidProjPoint ); - } - - return; -}