From 6bc94c22113c7133067ebd33f85ccb437b3c3fe0 Mon Sep 17 00:00:00 2001 From: eap Date: Fri, 16 Oct 2020 01:16:47 +0300 Subject: [PATCH] bos #20144 [CEA 20142] Import1D - Error: Algorithm failed --- src/SMESHUtils/SMESH_MAT2d.cxx | 5 + src/SMESHUtils/SMESH_MeshAlgos.cxx | 2 +- src/SMESHUtils/SMESH_Octree.cxx | 10 +- src/SMESHUtils/SMESH_Tree.hxx | 5 +- src/StdMeshers/StdMeshers_Hexa_3D.cxx | 10 - src/StdMeshers/StdMeshers_Import_1D.cxx | 389 +++++++++++++++++++--- src/StdMeshers/StdMeshers_Import_1D2D.cxx | 18 +- 7 files changed, 371 insertions(+), 68 deletions(-) diff --git a/src/SMESHUtils/SMESH_MAT2d.cxx b/src/SMESHUtils/SMESH_MAT2d.cxx index 601fdd84d..fbbb0553f 100644 --- a/src/SMESHUtils/SMESH_MAT2d.cxx +++ b/src/SMESHUtils/SMESH_MAT2d.cxx @@ -1077,6 +1077,11 @@ namespace } } + else // 2D_mesh_QuadranglePreference_00/A1, bos20144.brep + { + continue; // bndSegs.size() == 1 + } + bndSegs[i].setBranch( branchID, bndSegsPerEdge ); // set to i-th and to the opposite bndSeg if ( bndSegs[i].hasOppositeEdge() ) branchEdges[ bndSegs[i].branchID() ].push_back( bndSegs[i]._edge ); diff --git a/src/SMESHUtils/SMESH_MeshAlgos.cxx b/src/SMESHUtils/SMESH_MeshAlgos.cxx index 14c781836..db6f5e690 100644 --- a/src/SMESHUtils/SMESH_MeshAlgos.cxx +++ b/src/SMESHUtils/SMESH_MeshAlgos.cxx @@ -317,7 +317,7 @@ namespace // Utils used in SMESH_ElementSearcherImpl::FindElementsByPoint() //================================================================================ /*! - * \brief Redistrubute element boxes among children + * \brief Redistribute element boxes among children */ //================================================================================ diff --git a/src/SMESHUtils/SMESH_Octree.cxx b/src/SMESHUtils/SMESH_Octree.cxx index 19473e226..85df6d93b 100644 --- a/src/SMESHUtils/SMESH_Octree.cxx +++ b/src/SMESHUtils/SMESH_Octree.cxx @@ -28,6 +28,9 @@ // #include "SMESH_Octree.hxx" +#include +#include + //=========================================================================== /*! * Constructor. limit must be provided at tree root construction. @@ -90,7 +93,12 @@ void SMESH_Octree::enlargeByFactor( Bnd_B3d* box, double factor ) const { gp_XYZ halfSize = 0.5 * ( box->CornerMax() - box->CornerMin() ); for ( int iDim = 1; iDim <= 3; ++iDim ) - halfSize.SetCoord( iDim, factor * halfSize.Coord( iDim )); + { + double newHSize = factor * halfSize.Coord( iDim ); + if ( newHSize < std::numeric_limits::min() ) + newHSize = Precision::Confusion(); // 1.e-7 + halfSize.SetCoord( iDim, newHSize ); + } box->SetHSize( halfSize ); } } diff --git a/src/SMESHUtils/SMESH_Tree.hxx b/src/SMESHUtils/SMESH_Tree.hxx index 37e92a18c..7be283943 100644 --- a/src/SMESHUtils/SMESH_Tree.hxx +++ b/src/SMESHUtils/SMESH_Tree.hxx @@ -31,6 +31,8 @@ #include "SMESH_Utils.hxx" +const double theEnlargeFactor = 1. + 1e-10; + //================================================================================ // Data limiting the tree height struct SMESH_TreeLimit { @@ -160,6 +162,7 @@ void SMESH_Tree::compute() { if ( !myLimit ) myLimit = new SMESH_TreeLimit(); myBox = buildRootBox(); + enlargeByFactor( myBox, theEnlargeFactor ); if ( myLimit->myMinBoxSize > 0. && maxSize() <= myLimit->myMinBoxSize ) myIsLeaf = true; else @@ -227,7 +230,7 @@ void SMESH_Tree::buildChildren() myChildren[i]->myLimit = myLimit; myChildren[i]->myLevel = myLevel + 1; myChildren[i]->myBox = newChildBox( i ); - enlargeByFactor( myChildren[i]->myBox, 1. + 1e-10 ); + enlargeByFactor( myChildren[i]->myBox, theEnlargeFactor ); if ( myLimit->myMinBoxSize > 0. && myChildren[i]->maxSize() <= myLimit->myMinBoxSize ) myChildren[i]->myIsLeaf = true; } diff --git a/src/StdMeshers/StdMeshers_Hexa_3D.cxx b/src/StdMeshers/StdMeshers_Hexa_3D.cxx index cde18b142..db2f66880 100644 --- a/src/StdMeshers/StdMeshers_Hexa_3D.cxx +++ b/src/StdMeshers/StdMeshers_Hexa_3D.cxx @@ -840,16 +840,6 @@ bool StdMeshers_Hexa_3D::Compute(SMESH_Mesh & aMesh, vector< const SMDS_MeshNode* >& columnX0 = columns[ colIndex( x, 0 )]; for ( z = 0; z < zSize; ++z ) renumHelper.AddReplacingNode( columnX0[ z ] ); - if ( x == 0 || x == X ) - { - for ( y = 1; y < ySize; ++y ) - { - vector< const SMDS_MeshNode* >& column0Y = columns[ colIndex( x, y )]; - for ( z = 0; z < zSize; ++z ) - renumHelper.AddReplacingNode( column0Y[ z ] ); - } - continue; - } } const double rX = x / double(X); diff --git a/src/StdMeshers/StdMeshers_Import_1D.cxx b/src/StdMeshers/StdMeshers_Import_1D.cxx index 179dd52e4..e3d59a79c 100644 --- a/src/StdMeshers/StdMeshers_Import_1D.cxx +++ b/src/StdMeshers/StdMeshers_Import_1D.cxx @@ -38,14 +38,19 @@ #include "SMESH_Mesh.hxx" #include "SMESH_MeshEditor.hxx" #include "SMESH_MesherHelper.hxx" +#include "SMESH_Octree.hxx" #include "SMESH_subMesh.hxx" #include "SMESH_subMeshEventListener.hxx" #include "Utils_SALOME_Exception.hxx" #include "utilities.h" +#include #include #include +#include +#include +#include #include #include #include @@ -55,25 +60,251 @@ using namespace std; -//============================================================================= -/*! - * Creates StdMeshers_Import_1D - */ -//============================================================================= - -StdMeshers_Import_1D::StdMeshers_Import_1D(int hypId, SMESH_Gen * gen) - :SMESH_1D_Algo(hypId, gen), _sourceHyp(0) -{ - _name = "Import_1D"; - _shapeType = (1 << TopAbs_EDGE); - - _compatibleHypothesis.push_back("ImportSource1D"); -} - //================================================================================ namespace // INTERNAL STUFF //================================================================================ { + /*! + * \brief Compute point position on a curve. Use octree to fast reject far points + */ + class CurveProjector : public SMESH_Octree + { + public: + CurveProjector( const TopoDS_Edge& edge, double enlarge ); + + bool IsOnCurve( const gp_XYZ& point, double & distance2, double & u ); + + bool IsOut( const gp_XYZ& point ) const { return getBox()->IsOut( point ); } + + protected: + CurveProjector() {} + SMESH_Octree* newChild() const { return new CurveProjector; } + void buildChildrenData(); + Bnd_B3d* buildRootBox(); + + private: + struct CurveSegment : public Bnd_B3d + { + double _chord, _chord2, _length2; + gp_Pnt _pFirst, _pLast; + gp_Lin _line; + Handle(Geom_Curve) _curve; + + CurveSegment() {} + void Init( const gp_Pnt& pf, const gp_Pnt& pl, + double uf, double ul, double tol, Handle(Geom_Curve)& curve ); + bool IsOn( const gp_XYZ& point, double & distance2, double & u ); + bool IsInContact( const Bnd_B3d& bb ); + }; + std::vector< CurveSegment > _segments; + }; + + //=============================================================================== + /*! + * \brief Create an octree of curve segments + */ + //================================================================================ + + CurveProjector::CurveProjector( const TopoDS_Edge& edge, double enlarge ) + :SMESH_Octree( 0 ) + { + double f,l; + Handle(Geom_Curve) curve = BRep_Tool::Curve( edge, f, l ); + double curDeflect = 0.3; // Curvature deflection + double angDeflect = 1e+100; // Angular deflection - don't control chordal error + GCPnts_TangentialDeflection div( BRepAdaptor_Curve( edge ), angDeflect, curDeflect ); + _segments.resize( div.NbPoints() - 1 ); + for ( int i = 1; i < div.NbPoints(); ++i ) + try { + _segments[ i - 1 ].Init( div.Value( i ), div.Value( i+1 ), + div.Parameter( i ), div.Parameter( i+1 ), + enlarge, curve ); + } + catch ( Standard_Failure ) { + _segments.resize( _segments.size() - 1 ); + --i; + } + if ( _segments.size() < 3 ) + myIsLeaf = true; + + compute(); + + if ( _segments.size() == 1 ) + myBox->Enlarge( enlarge ); + } + + //================================================================================ + /*! + * \brief Return the maximal box + */ + //================================================================================ + + Bnd_B3d* CurveProjector::buildRootBox() + { + Bnd_B3d* box = new Bnd_B3d; + for ( size_t i = 0; i < _segments.size(); ++i ) + box->Add( _segments[i] ); + return box; + } + + //================================================================================ + /*! + * \brief Redistribute segments among children + */ + //================================================================================ + + void CurveProjector::buildChildrenData() + { + bool allIn = true; + for ( size_t i = 0; i < _segments.size(); ++i ) + { + for (int j = 0; j < 8; j++) + { + if ( _segments[i].IsInContact( *myChildren[j]->getBox() )) + ((CurveProjector*)myChildren[j])->_segments.push_back( _segments[i]); + else + allIn = false; + } + } + if ( allIn && _segments.size() < 3 ) + { + myIsLeaf = true; + for (int j = 0; j < 8; j++) + static_cast( myChildren[j])->myIsLeaf = true; + } + else + { + SMESHUtils::FreeVector( _segments ); // = _segments.clear() + free memory + + for (int j = 0; j < 8; j++) + { + CurveProjector* child = static_cast( myChildren[j]); + if ( child->_segments.size() < 3 ) + child->myIsLeaf = true; + } + } + } + + //================================================================================ + /*! + * \brief Return true if a point is close to the curve + * \param [in] point - the point + * \param [out] distance2 - distance to the curve + * \param [out] u - parameter on the curve + * \return bool - is the point is close to the curve + */ + //================================================================================ + + bool CurveProjector::IsOnCurve( const gp_XYZ& point, double & distance2, double & u ) + { + if ( getBox()->IsOut( point )) + return false; + + if ( isLeaf() ) + { + for ( size_t i = 0; i < _segments.size(); ++i ) + if ( !_segments[i].IsOut( point ) && + _segments[i].IsOn( point, distance2, u )) + return true; + } + else + { + for (int i = 0; i < 8; i++) + if (((CurveProjector*) myChildren[i])->IsOnCurve( point, distance2, u )) + return true; + } + return false; + } + + //================================================================================ + /*! + * \brief Initialize + */ + //================================================================================ + + void CurveProjector::CurveSegment::Init(const gp_Pnt& pf, + const gp_Pnt& pl, + const double uf, + const double ul, + const double tol, + Handle(Geom_Curve)& curve ) + { + _pFirst = pf; + _pLast = pl; + _curve = curve; + _length2 = pf.SquareDistance( pl ); + _chord2 = Max( _line. SquareDistance( curve->Value( uf + 0.25 * ( ul - uf ))), + Max( _line.SquareDistance( curve->Value( uf + 0.5 * ( ul - uf ))), + _line.SquareDistance( curve->Value( uf + 0.75 * ( ul - uf ))))); + _chord2 = Max( tol, _chord2 ); + _chord = Sqrt( _chord2 ); + _line.SetLocation( pf ); + _line.SetDirection( gp_Vec( pf, pl )); + + Bnd_Box bb; + BndLib_Add3dCurve::Add( GeomAdaptor_Curve( curve, uf, ul ), tol, bb ); + Add( bb.CornerMin() ); + Add( bb.CornerMax() ); + } + + //================================================================================ + /*! + * \brief Return true if a point is close to the curve segment + * \param [in] point - the point + * \param [out] distance2 - distance to the curve + * \param [out] u - parameter on the curve + * \return bool - is the point is close to the curve segment + */ + //================================================================================ + + bool CurveProjector::CurveSegment::IsOn( const gp_XYZ& point, double & distance2, double & u ) + { + distance2 = _line.SquareDistance( point ); + if ( distance2 > _chord2 ) + return false; + + // check if the point projection falls into the segment range + { + gp_Vec edge( _pFirst, _pLast ); + gp_Vec n1p ( _pFirst, point ); + u = ( edge * n1p ) / _length2; // param [0,1] on the edge + if ( u < 0 ) + { + if ( _pFirst.SquareDistance( point ) > _chord2 ) + return false; + } + else if ( u > _chord ) + { + if ( _pLast.SquareDistance( point ) > _chord2 ) + return false; + } + } + gp_Pnt proj; + distance2 = ShapeAnalysis_Curve().Project( _curve, point, Precision::Confusion(), + proj, u, false ); + distance2 *= distance2; + return true; + } + + //================================================================================ + /*! + * \brief Check if the segment is in contact with a box + */ + //================================================================================ + + bool CurveProjector::CurveSegment::IsInContact( const Bnd_B3d& bb ) + { + if ( bb.IsOut( _line.Position(), /*isRay=*/true, _chord )) + return false; + + gp_Ax1 axRev = _line.Position().Reversed(); + axRev.SetLocation( _pLast ); + return !bb.IsOut( axRev, /*isRay=*/true, _chord ); + } + + //================================================================================ + //================================================================================ + int getSubmeshIDForCopiedMesh(const SMESHDS_Mesh* srcMeshDS, SMESH_Mesh* tgtMesh); enum _ListenerDataType @@ -590,8 +821,50 @@ namespace // INTERNAL STUFF return 0; } + //================================================================================ + /*! + * \brief Return minimal square length of edges of 1D and 2D elements sharing the node + */ + //================================================================================ + + double getMinEdgeLength2( const SMDS_MeshNode* n ) + { + SMESH_NodeXYZ p = n; + double minLen2 = Precision::Infinite(); + for ( SMDS_ElemIteratorPtr eIt = n->GetInverseElementIterator(); eIt->more(); ) + { + const SMDS_MeshElement* e = eIt->next(); + const SMDSAbs_ElementType type = e->GetType(); + if ( type != SMDSAbs_Edge && type != SMDSAbs_Face ) + continue; + int i = e->GetNodeIndex( n ); + int iNext = SMESH_MesherHelper::WrapIndex( i + 1, e->NbCornerNodes() ); + minLen2 = Min( minLen2, p.SquareDistance( e->GetNode( iNext ))); + if ( type != SMDSAbs_Face ) + continue; + int iPrev = SMESH_MesherHelper::WrapIndex( i - 1, e->NbCornerNodes() ); + minLen2 = Min( minLen2, p.SquareDistance( e->GetNode( iPrev ))); + } + return minLen2; + } + } // namespace +//============================================================================= +/*! + * Creates StdMeshers_Import_1D + */ +//============================================================================= + +StdMeshers_Import_1D::StdMeshers_Import_1D(int hypId, SMESH_Gen * gen) + :SMESH_1D_Algo(hypId, gen), _sourceHyp(0) +{ + _name = "Import_1D"; + _shapeType = (1 << TopAbs_EDGE); + + _compatibleHypothesis.push_back("ImportSource1D"); +} + //============================================================================= /*! * Check presence of a hypothesis @@ -658,17 +931,31 @@ bool StdMeshers_Import_1D::Compute(SMESH_Mesh & theMesh, const TopoDS_Shape & th const double edgeTol = BRep_Tool::Tolerance( geomEdge ); const int shapeID = tgtMesh->ShapeToIndex( geomEdge ); - set subShapeIDs; - subShapeIDs.insert( shapeID ); + + double geomTol = Precision::Confusion(); + for ( size_t iG = 0; iG < srcGroups.size(); ++iG ) + { + const SMESHDS_GroupBase* srcGroup = srcGroups[iG]->GetGroupDS(); + for ( SMDS_ElemIteratorPtr srcElems = srcGroup->GetElements(); srcElems->more(); ) + { + const SMDS_MeshElement* edge = srcElems->next(); + geomTol = Sqrt( 0.5 * ( getMinEdgeLength2( edge->GetNode(0) ) + + getMinEdgeLength2( edge->GetNode(1) ))) / 25; + iG = srcGroups.size(); + break; + } + } + CurveProjector curveProjector( geomEdge, geomTol ); // get nodes on vertices + set vertexIDs; list < SMESH_TNodeXYZ > vertexNodes; list < SMESH_TNodeXYZ >::iterator vNIt; TopExp_Explorer vExp( theShape, TopAbs_VERTEX ); for ( ; vExp.More(); vExp.Next() ) { const TopoDS_Vertex& v = TopoDS::Vertex( vExp.Current() ); - if ( !subShapeIDs.insert( tgtMesh->ShapeToIndex( v )).second ) + if ( !vertexIDs.insert( tgtMesh->ShapeToIndex( v )).second ) continue; // closed edge const SMDS_MeshNode* n = SMESH_Algo::VertexNode( v, tgtMesh ); if ( !n ) @@ -696,56 +983,61 @@ bool StdMeshers_Import_1D::Compute(SMESH_Mesh & theMesh, const TopoDS_Shape & th SMDS_ElemIteratorPtr srcElems = srcGroup->GetElements(); vector newNodes; - SMDS_MeshNode *tmpNode = helper.AddNode(0,0,0); - double u = 0.314159; // "random" value between 0 and 1, avoid 0 and 1, false detection possible on edge restrictions while ( srcElems->more() ) // loop on group contents { const SMDS_MeshElement* edge = srcElems->next(); + gp_XYZ middle = 0.5 * ( SMESH_NodeXYZ( edge->GetNode(0)) + + SMESH_NodeXYZ( edge->GetNode(1))); + if ( curveProjector.IsOut( middle )) + continue; + // find or create nodes of a new edge newNodes.resize( edge->NbNodes() ); - //MESSAGE("edge->NbNodes " << edge->NbNodes()); newNodes.back() = 0; + int nbNodesOnVertex = 0; SMDS_MeshElement::iterator node = edge->begin_nodes(); - SMESH_TNodeXYZ a(edge->GetNode(0)); - // --- define a tolerance relative to the length of an edge - double mytol = a.Distance(edge->GetNode(edge->NbNodes()-1))/25; - //mytol = max(1.E-5, 10*edgeTol); // too strict and not necessary - //MESSAGE("mytol = " << mytol); for ( size_t i = 0; i < newNodes.size(); ++i, ++node ) { - TNodeNodeMap::iterator n2nIt = n2n->insert( make_pair( *node, (SMDS_MeshNode*)0 )).first; + TNodeNodeMap::iterator n2nIt = n2n->insert( make_pair( *node, nullptr )).first; if ( n2nIt->second ) { - if ( !subShapeIDs.count( n2nIt->second->getshapeId() )) - break; + int sId = n2nIt->second->getshapeId(); + if ( sId != shapeID ) + { + if ( vertexIDs.count( sId )) + ++nbNodesOnVertex; + else + break; + } } - else + else if ( !vertexNodes.empty() ) { // find an existing vertex node double checktol = max(1.E-10, 10*edgeTol*edgeTol); for ( vNIt = vertexNodes.begin(); vNIt != vertexNodes.end(); ++vNIt) if ( vNIt->SquareDistance( *node ) < checktol) { - //MESSAGE("SquareDistance " << vNIt->SquareDistance( *node ) << " checktol " << checktol <<" "<X()<<" "<Y()<<" "<Z()); (*n2nIt).second = vNIt->_node; vertexNodes.erase( vNIt ); + ++nbNodesOnVertex; break; } - else if ( vNIt->SquareDistance( *node ) < 10*checktol) - MESSAGE("SquareDistance missed" << vNIt->SquareDistance( *node ) << " checktol " << checktol <<" "<X()<<" "<Y()<<" "<Z()); } if ( !n2nIt->second ) { - // find out if node lies on theShape - //double dxyz[4]; - tmpNode->setXYZ( (*node)->X(), (*node)->Y(), (*node)->Z()); - if ( helper.CheckNodeU( geomEdge, tmpNode, u, mytol, /*force=*/true)) // , dxyz )) // dxyz used for debug purposes + // find out if the node lies on theShape + SMESH_NodeXYZ xyz = *node; + double dist2, u; + if ( curveProjector.IsOnCurve( xyz, dist2, u )) { - SMDS_MeshNode* newNode = tgtMesh->AddNode( (*node)->X(), (*node)->Y(), (*node)->Z()); - n2nIt->second = newNode; - tgtMesh->SetNodeOnEdge( newNode, shapeID, u ); - //MESSAGE("u=" << u << " " << newNode->X()<< " " << newNode->Y()<< " " << newNode->Z()); - //MESSAGE("d=" << dxyz[0] << " " << dxyz[1] << " " << dxyz[2] << " " << dxyz[3]); + // tolerance relative to the length of surrounding edges + double mytol2 = getMinEdgeLength2( *node ) / 25 / 25; + if ( dist2 < mytol2 ) + { + SMDS_MeshNode* newNode = tgtMesh->AddNode( xyz.X(), xyz.Y(), xyz.Z() ); + n2nIt->second = newNode; + tgtMesh->SetNodeOnEdge( newNode, shapeID, u ); + } } } if ( !(newNodes[i] = n2nIt->second )) @@ -763,12 +1055,17 @@ bool StdMeshers_Import_1D::Compute(SMESH_Mesh & theMesh, const TopoDS_Shape & th newEdge = tgtMesh->AddEdge( newNodes[0], newNodes[1], newNodes[2] ); else newEdge = tgtMesh->AddEdge( newNodes[0], newNodes[1]); - //MESSAGE("add Edge " << newNodes[0]->GetID() << " " << newNodes[1]->GetID()); tgtMesh->SetMeshElementOnShape( newEdge, shapeID ); e2e->insert( make_pair( edge, newEdge )); - } - helper.GetMeshDS()->RemoveNode(tmpNode); - } + + if ( nbNodesOnVertex >= 2 ) // EDGE is meshed by a sole segment + { + iG = srcGroups.size(); // stop looingp on groups + break; + } + } // loop on group contents + } // loop on groups + if ( n2n->empty()) return error("Empty source groups"); diff --git a/src/StdMeshers/StdMeshers_Import_1D2D.cxx b/src/StdMeshers/StdMeshers_Import_1D2D.cxx index 04d903508..30f3d163c 100644 --- a/src/StdMeshers/StdMeshers_Import_1D2D.cxx +++ b/src/StdMeshers/StdMeshers_Import_1D2D.cxx @@ -324,8 +324,7 @@ bool StdMeshers_Import_1D2D::Compute(SMESH_Mesh & theMesh, const TopoDS_Shape & { const SMDS_MeshElement* face = srcElems->next(); - SMDS_MeshElement::iterator node = face->begin_nodes(); - if ( bndBox3d.IsOut( SMESH_TNodeXYZ( *node ))) + if ( bndBox3d.IsOut( SMESH_NodeXYZ( face->GetNode(0) ))) continue; // find or create nodes of a new face @@ -334,13 +333,14 @@ bool StdMeshers_Import_1D2D::Compute(SMESH_Mesh & theMesh, const TopoDS_Shape & newNodes.back() = 0; int nbCreatedNodes = 0; bool isOut = false, isIn = false; // if at least one node isIn - do not classify other nodes - for ( size_t i = 0; i < newNodes.size(); ++i, ++node ) + for ( size_t i = 0; i < newNodes.size(); ++i ) { - SMESH_TNodeXYZ nXYZ = *node; + const SMDS_MeshNode* node = face->GetNode( i ); + SMESH_NodeXYZ nXYZ = node; nodeState[ i ] = TopAbs_UNKNOWN; newNodes [ i ] = 0; - it_isnew = n2n->insert( make_pair( *node, (SMDS_MeshNode*)0 )); + it_isnew = n2n->insert( make_pair( node, nullptr )); n2nIt = it_isnew.first; const SMDS_MeshNode* & newNode = n2nIt->second; @@ -354,7 +354,7 @@ bool StdMeshers_Import_1D2D::Compute(SMESH_Mesh & theMesh, const TopoDS_Shape & if ( newNode->GetID() < (int) isNodeIn.size() && isNodeIn[ newNode->GetID() ]) isIn = true; - if ( !isIn && bndNodes.count( *node )) + if ( !isIn && bndNodes.count( node )) nodeState[ i ] = TopAbs_ON; } else @@ -373,7 +373,7 @@ bool StdMeshers_Import_1D2D::Compute(SMESH_Mesh & theMesh, const TopoDS_Shape & { // find out if node lies on the surface of theShape gp_XY uv( Precision::Infinite(), 0 ); - isOut = ( !helper.CheckNodeUV( geomFace, *node, uv, groupTol, /*force=*/true ) || + isOut = ( !helper.CheckNodeUV( geomFace, node, uv, groupTol, /*force=*/true ) || bndBox2d.IsOut( uv )); //int iCoo; if ( !isOut && !isIn ) // classify @@ -381,7 +381,7 @@ bool StdMeshers_Import_1D2D::Compute(SMESH_Mesh & theMesh, const TopoDS_Shape & nodeState[i] = classifier.Perform( uv ); //classifier.Perform( geomFace, uv, clsfTol ); //nodeState[i] = classifier.State(); isOut = ( nodeState[i] == TopAbs_OUT ); - if ( isOut && helper.IsOnSeam( uv ) && onEdgeClassifier.IsSatisfy( (*node)->GetID() )) + if ( isOut && helper.IsOnSeam( uv ) && onEdgeClassifier.IsSatisfy( node->GetID() )) { // uv.SetCoord( iCoo, helper.GetOtherParam( uv.Coord( iCoo ))); // classifier.Perform( geomFace, uv, clsfTol ); @@ -402,7 +402,7 @@ bool StdMeshers_Import_1D2D::Compute(SMESH_Mesh & theMesh, const TopoDS_Shape & isNodeIn.resize( newNode->GetID() + 1, false ); } if ( nodeState[i] == TopAbs_ON ) - bndNodes.insert( *node ); + bndNodes.insert( node ); else if ( nodeState[i] != TopAbs_UNKNOWN ) isNodeIn[ newNode->GetID() ] = isIn = true; } -- 2.30.2