From 4466dfe1ce1644260fb11df461175e53785fb867 Mon Sep 17 00:00:00 2001 From: eap Date: Tue, 18 Jan 2011 12:09:32 +0000 Subject: [PATCH] 0020832: EDF 1359 SMESH : Automatic meshing of boundary layers 1) make CheckNodeUV() and CheckNodeU() optionally return XYZ of node projection to shape 2) prevent ancestors iterator from returning duplicates 3) move IsClosedEdge() to SMESH_MesherHelper from StdMeshers_ProjectionUtils 4) for hexa 3D static bool LoadNodeColumns(TParam2ColumnMap & theParam2ColumnMap, const TopoDS_Face& theFace, const TopoDS_Edge& theBaseEdge, SMESHDS_Mesh* theMesh, + SMESH_ProxyMesh* theProxyMesh=0); --- src/SMESH/SMESH_MesherHelper.cxx | 128 +++++++++++++++++++++++-------- src/SMESH/SMESH_MesherHelper.hxx | 15 +++- 2 files changed, 109 insertions(+), 34 deletions(-) diff --git a/src/SMESH/SMESH_MesherHelper.cxx b/src/SMESH/SMESH_MesherHelper.cxx index 3c3433548..66e630f0a 100644 --- a/src/SMESH/SMESH_MesherHelper.cxx +++ b/src/SMESH/SMESH_MesherHelper.cxx @@ -30,6 +30,7 @@ #include "SMDS_EdgePosition.hxx" #include "SMDS_VolumeTool.hxx" #include "SMESH_subMesh.hxx" +#include "SMESH_ProxyMesh.hxx" #include #include @@ -291,7 +292,7 @@ bool SMESH_MesherHelper::IsMedium(const SMDS_MeshNode* node, //======================================================================= TopoDS_Shape SMESH_MesherHelper::GetSubShapeByNode(const SMDS_MeshNode* node, - SMESHDS_Mesh* meshDS) + const SMESHDS_Mesh* meshDS) { int shapeID = node->getshapeId(); if ( 0 < shapeID && shapeID <= meshDS->MaxShapeIndex() ) @@ -500,7 +501,8 @@ bool SMESH_MesherHelper::CheckNodeUV(const TopoDS_Face& F, const SMDS_MeshNode* n, gp_XY& uv, const double tol, - const bool force) const + const bool force, + double distXYZ[4]) const { int shapeID = n->getshapeId(); if ( force || toCheckPosOnShape( shapeID )) @@ -511,13 +513,19 @@ bool SMESH_MesherHelper::CheckNodeUV(const TopoDS_Face& F, // check that uv is correct TopLoc_Location loc; Handle(Geom_Surface) surface = BRep_Tool::Surface( F,loc ); - gp_Pnt nodePnt = XYZ( n ); + gp_Pnt nodePnt = XYZ( n ), surfPnt(0,0,0); + double dist = 0; if ( !loc.IsIdentity() ) nodePnt.Transform( loc.Transformation().Inverted() ); - if ( Precision::IsInfinite( uv.X() ) || - Precision::IsInfinite( uv.Y() ) || - nodePnt.Distance( surface->Value( uv.X(), uv.Y() )) > toldis ) + bool infinit = ( Precision::IsInfinite( uv.X() ) || Precision::IsInfinite( uv.Y() )); + if ( infinit || + (dist = nodePnt.Distance( surfPnt = surface->Value( uv.X(), uv.Y() ))) > tol ) { setPosOnShapeValidity( shapeID, false ); + if ( !infinit && distXYZ ) { + surfPnt.Transform( loc ); + distXYZ[0] = dist; + distXYZ[1] = surfPnt.X(); distXYZ[2] = surfPnt.Y(); distXYZ[3]=surfPnt.Z(); + } // uv incorrect, project the node to surface GeomAPI_ProjectPointOnSurf& projector = GetProjector( F, loc, tol ); projector.Perform( nodePnt ); @@ -529,7 +537,14 @@ bool SMESH_MesherHelper::CheckNodeUV(const TopoDS_Face& F, Quantity_Parameter U,V; projector.LowerDistanceParameters(U,V); uv.SetCoord( U,V ); - if ( nodePnt.Distance( surface->Value( U, V )) > toldis ) + surfPnt = surface->Value( U, V ); + dist = nodePnt.Distance( surfPnt ); + if ( distXYZ ) { + surfPnt.Transform( loc ); + distXYZ[0] = dist; + distXYZ[1] = surfPnt.X(); distXYZ[2] = surfPnt.Y(); distXYZ[3]=surfPnt.Z(); + } + if ( dist > tol ) { MESSAGE( "SMESH_MesherHelper::CheckNodeUV(), invalid projection" ); return false; @@ -689,14 +704,14 @@ bool SMESH_MesherHelper::CheckNodeU(const TopoDS_Edge& E, double& u, const double tol, const bool force, - double* distance) const + double distXYZ[4]) const { int shapeID = n->getshapeId(); if ( force || toCheckPosOnShape( shapeID )) { - double toldis = tol; - double tolmin = 1.e-7*myMesh->GetMeshDS()->getMaxDim(); // nodes coordinates are stored in float format - if (toldis < tolmin) toldis = tolmin; + //double toldis = tol; + //double tolmin = 1.e-7*myMesh->GetMeshDS()->getMaxDim(); // nodes coordinates are stored in float format + //if (toldis < tolmin) toldis = tolmin; // check that u is correct TopLoc_Location loc; double f,l; Handle(Geom_Curve) curve = BRep_Tool::Curve( E,loc,f,l ); @@ -712,9 +727,14 @@ bool SMESH_MesherHelper::CheckNodeU(const TopoDS_Edge& E, { gp_Pnt nodePnt = SMESH_MeshEditor::TNodeXYZ( n ); if ( !loc.IsIdentity() ) nodePnt.Transform( loc.Transformation().Inverted() ); - double dist = nodePnt.Distance( curve->Value( u )); - if ( distance ) *distance = dist; - if ( dist > toldis ) + gp_Pnt curvPnt = curve->Value( u ); + double dist = nodePnt.Distance( curvPnt ); + if ( distXYZ ) { + curvPnt.Transform( loc ); + distXYZ[0] = dist; + distXYZ[1] = curvPnt.X(); distXYZ[2] = curvPnt.Y(); distXYZ[3]=curvPnt.Z(); + } + if ( dist > tol /*toldis*/ ) { setPosOnShapeValidity( shapeID, false ); // u incorrect, project the node to the curve @@ -736,12 +756,17 @@ bool SMESH_MesherHelper::CheckNodeU(const TopoDS_Edge& E, } Quantity_Parameter U = projector->LowerDistanceParameter(); u = double( U ); - dist = nodePnt.Distance( curve->Value( U )); - if ( distance ) *distance = dist; - if ( dist > toldis ) + curvPnt = curve->Value( u ); + dist = nodePnt.Distance( curvPnt ); + if ( distXYZ ) { + curvPnt.Transform( loc ); + distXYZ[0] = dist; + distXYZ[1] = curvPnt.X(); distXYZ[2] = curvPnt.Y(); distXYZ[3]=curvPnt.Z(); + } + if ( dist > tol /*toldis*/) { MESSAGE( "SMESH_MesherHelper::CheckNodeU(), invalid projection" ); - MESSAGE("distance " << nodePnt.Distance(curve->Value( U )) << " " << toldis); + MESSAGE("distance " << dist << " " << tol/*dis*/); return false; } // store the fixed U on the edge @@ -948,7 +973,7 @@ const SMDS_MeshNode* SMESH_MesherHelper::getMediumNodeOnComposedWire(const SMDS_ // To find position on edge and 3D position for n12, // project to 2 edges and select projection most close to - double u = 0, distMiddleProj = Precision::Infinite(); + double u = 0, distMiddleProj = Precision::Infinite(), distXYZ[4]; int iOkEdge = 0; TopoDS_Edge edges[2]; for ( int is2nd = 0; is2nd < 2; ++is2nd ) @@ -962,11 +987,11 @@ const SMDS_MeshNode* SMESH_MesherHelper::getMediumNodeOnComposedWire(const SMDS_ // project to get U of projection and distance from middle to projection TopoDS_Edge edge = edges[ is2nd ] = TopoDS::Edge( shape ); double node2MiddleDist = middle.Distance( XYZ(n) ); - double foundU = GetNodeU( edge, n ), foundDist = node2MiddleDist; - CheckNodeU( edge, n12, foundU, 2*BRep_Tool::Tolerance(edge), /*force=*/true, &foundDist ); - if ( foundDist < node2MiddleDist ) + double foundU = GetNodeU( edge, n ); + CheckNodeU( edge, n12, foundU, 2*BRep_Tool::Tolerance(edge), /*force=*/true, distXYZ ); + if ( distXYZ[0] < node2MiddleDist ) { - distMiddleProj = foundDist; + distMiddleProj = distXYZ[0]; u = foundU; iOkEdge = is2nd; } @@ -1452,9 +1477,24 @@ SMESH_MesherHelper::AddPolyhedralVolume (const std::vector bool SMESH_MesherHelper::LoadNodeColumns(TParam2ColumnMap & theParam2ColumnMap, const TopoDS_Face& theFace, const TopoDS_Edge& theBaseEdge, - SMESHDS_Mesh* theMesh) + SMESHDS_Mesh* theMesh, + SMESH_ProxyMesh* theProxyMesh) { - SMESHDS_SubMesh* faceSubMesh = theMesh->MeshElements( theFace ); + const SMESHDS_SubMesh* faceSubMesh = 0; + if ( theProxyMesh ) + { + faceSubMesh = theProxyMesh->GetSubMesh( theFace ); + if ( !faceSubMesh || + faceSubMesh->NbElements() == 0 || + theProxyMesh->IsTemporary( faceSubMesh->GetElements()->next() )) + { + // can use a proxy sub-mesh with not temporary elements only + faceSubMesh = 0; + theProxyMesh = 0; + } + } + if ( !faceSubMesh ) + faceSubMesh = theMesh->MeshElements( theFace ); if ( !faceSubMesh || faceSubMesh->NbElements() == 0 ) return false; @@ -1475,12 +1515,21 @@ bool SMESH_MesherHelper::LoadNodeColumns(TParam2ColumnMap & theParam2ColumnMap, nCol.resize( nbRows ); nCol[0] = u_n->second; } + TParam2ColumnMap::iterator par_nVec_2, par_nVec_1 = theParam2ColumnMap.begin(); + if ( theProxyMesh ) + { + for ( ; par_nVec_1 != theParam2ColumnMap.end(); ++par_nVec_1 ) + { + const SMDS_MeshNode* & n = par_nVec_1->second[0]; + n = theProxyMesh->GetProxyNode( n ); + } + } // fill theParam2ColumnMap column by column by passing from nodes on // theBaseEdge up via mesh faces on theFace - TParam2ColumnMap::iterator par_nVec_2 = theParam2ColumnMap.begin(); - TParam2ColumnMap::iterator par_nVec_1 = par_nVec_2++; + par_nVec_2 = theParam2ColumnMap.begin(); + par_nVec_1 = par_nVec_2++; TIDSortedElemSet emptySet, avoidSet; for ( ; par_nVec_2 != theParam2ColumnMap.end(); ++par_nVec_1, ++par_nVec_2 ) { @@ -1511,7 +1560,7 @@ bool SMESH_MesherHelper::LoadNodeColumns(TParam2ColumnMap & theParam2ColumnMap, if ( iRow + 1 < nbRows ) // compact if necessary nCol1.resize( iRow + 1 ), nCol2.resize( iRow + 1 ); } - return true; + return theParam2ColumnMap.size() > 1 && theParam2ColumnMap.begin()->second.size() > 1; } //======================================================================= @@ -1611,6 +1660,21 @@ double SMESH_MesherHelper::MaxTolerance( const TopoDS_Shape& shape ) return tol; } +//================================================================================ +/*! + * \brief Check if the first and last vertices of an edge are the same + * \param anEdge - the edge to check + * \retval bool - true if same + */ +//================================================================================ + +bool SMESH_MesherHelper::IsClosedEdge( const TopoDS_Edge& anEdge ) +{ + if ( anEdge.Orientation() >= TopAbs_INTERNAL ) + return IsClosedEdge( TopoDS::Edge( anEdge.Oriented( TopAbs_FORWARD ))); + return TopExp::FirstVertex( anEdge ).IsSame( TopExp::LastVertex( anEdge )); +} + //======================================================================= //function : IsQuadraticMesh //purpose : Check mesh without geometry for: if all elements on this shape are quadratic, @@ -2932,10 +2996,14 @@ struct TAncestorsIterator : public SMDS_Iterator { TopTools_ListIteratorOfListOfShape _ancIter; TopAbs_ShapeEnum _type; + TopTools_MapOfShape _encountered; TAncestorsIterator( const TopTools_ListOfShape& ancestors, TopAbs_ShapeEnum type) : _ancIter( ancestors ), _type( type ) { - if ( _ancIter.More() && _ancIter.Value().ShapeType() != _type ) next(); + if ( _ancIter.More() ) { + if ( _ancIter.Value().ShapeType() != _type ) next(); + else _encountered.Add( _ancIter.Value() ); + } } virtual bool more() { @@ -2946,7 +3014,7 @@ struct TAncestorsIterator : public SMDS_Iterator const TopoDS_Shape* s = _ancIter.More() ? & _ancIter.Value() : 0; if ( _ancIter.More() ) for ( _ancIter.Next(); _ancIter.More(); _ancIter.Next()) - if ( _ancIter.Value().ShapeType() == _type ) + if ( _ancIter.Value().ShapeType() == _type && _encountered.Add( _ancIter.Value() )) break; return s; } diff --git a/src/SMESH/SMESH_MesherHelper.hxx b/src/SMESH/SMESH_MesherHelper.hxx index d904f95de..7a53c71eb 100644 --- a/src/SMESH/SMESH_MesherHelper.hxx +++ b/src/SMESH/SMESH_MesherHelper.hxx @@ -43,6 +43,7 @@ class GeomAPI_ProjectPointOnSurf; class GeomAPI_ProjectPointOnCurve; +class SMESH_ProxyMesh; typedef std::map TLinkNodeMap; typedef std::map::iterator ItTLinkNode; @@ -96,7 +97,8 @@ public: static bool LoadNodeColumns(TParam2ColumnMap & theParam2ColumnMap, const TopoDS_Face& theFace, const TopoDS_Edge& theBaseEdge, - SMESHDS_Mesh* theMesh); + SMESHDS_Mesh* theMesh, + SMESH_ProxyMesh* theProxyMesh=0); /*! * \brief Return support shape of a node * \param node - the node @@ -104,7 +106,7 @@ public: * \retval TopoDS_Shape - found support shape */ static TopoDS_Shape GetSubShapeByNode(const SMDS_MeshNode* node, - SMESHDS_Mesh* meshDS); + const SMESHDS_Mesh* meshDS); /*! * \brief Return a valid node index, fixing the given one if necessary @@ -143,6 +145,8 @@ public: static double MaxTolerance( const TopoDS_Shape& shape ); + static bool IsClosedEdge( const TopoDS_Edge& anEdge ); + public: // ---------- PUBLIC INSTANCE METHODS ---------- @@ -301,16 +305,19 @@ public: /*! * \brief Check and fix node UV on a face * \param force - check even if checks of other nodes on this face passed OK + * \param distXYZ - returns result distance and point coordinates * \retval bool - false if UV is bad and could not be fixed */ bool CheckNodeUV(const TopoDS_Face& F, const SMDS_MeshNode* n, gp_XY& uv, const double tol, - const bool force=false) const; + const bool force=false, + double distXYZ[4]=0) const; /*! * \brief Check and fix node U on an edge * \param force - check even if checks of other nodes on this edge passed OK + * \param distXYZ - returns result distance and point coordinates * \retval bool - false if U is bad and could not be fixed */ bool CheckNodeU(const TopoDS_Edge& E, @@ -318,7 +325,7 @@ public: double& u, const double tol, const bool force=false, - double* distance=0) const; + double distXYZ[4]=0) const; /*! * \brief Return middle UV taking in account surface period */ -- 2.39.2