X-Git-Url: http://git.salome-platform.org/gitweb/?a=blobdiff_plain;f=src%2FStdMeshers%2FStdMeshers_Cartesian_3D.cxx;h=dc82aa09a035697663a0f87732b4048c945036b1;hb=63879abd44911e419cdfe602c837c8884f3fbe12;hp=ea641c2b118c0e8f021d18ceed0a065ae4e5c435;hpb=650885c1351fd0ec14226de37a67dec62329d83f;p=modules%2Fsmesh.git diff --git a/src/StdMeshers/StdMeshers_Cartesian_3D.cxx b/src/StdMeshers/StdMeshers_Cartesian_3D.cxx index ea641c2b1..dc82aa09a 100644 --- a/src/StdMeshers/StdMeshers_Cartesian_3D.cxx +++ b/src/StdMeshers/StdMeshers_Cartesian_3D.cxx @@ -1,4 +1,4 @@ -// Copyright (C) 2007-2019 CEA/DEN, EDF R&D, OPEN CASCADE +// Copyright (C) 2007-2021 CEA/DEN, EDF R&D, OPEN CASCADE // // Copyright (C) 2003-2007 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN, // CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS @@ -79,6 +79,7 @@ #include #include #include +#include #include #include #include @@ -97,7 +98,11 @@ #include -//#undef WITH_TBB +#ifdef _DEBUG_ +// #define _MY_DEBUG_ +// #undef WITH_TBB +#endif + #ifdef WITH_TBB #ifdef WIN32 @@ -114,10 +119,6 @@ using namespace std; using namespace SMESH; -#ifdef _DEBUG_ -//#define _MY_DEBUG_ -#endif - //============================================================================= /*! * Constructor @@ -171,6 +172,8 @@ namespace { typedef int TGeomID; // IDs of sub-shapes + const TGeomID theUndefID = 1e+9; + //============================================================================= // Definitions of internal utils // -------------------------------------------------------------------------- @@ -182,23 +185,47 @@ namespace Trans_INTERNAL // for INTERNAL FACE }; // -------------------------------------------------------------------------- + /*! + * \brief Sub-entities of a FACE neighboring its concave VERTEX. + * Help to avoid linking nodes on EDGEs that seem connected + * by the concave FACE but the link actually lies outside the FACE + */ + struct ConcaveFace + { + TGeomID _concaveFace; + TGeomID _edge1, _edge2; + TGeomID _v1, _v2; + ConcaveFace( int f=0, int e1=0, int e2=0, int v1=0, int v2=0 ) + : _concaveFace(f), _edge1(e1), _edge2(e2), _v1(v1), _v2(v2) {} + bool HasEdge( TGeomID edge ) const { return edge == _edge1 || edge == _edge2; } + bool HasVertex( TGeomID v ) const { return v == _v1 || v == _v2; } + void SetEdge( TGeomID edge ) { ( _edge1 ? _edge2 : _edge1 ) = edge; } + void SetVertex( TGeomID v ) { ( _v1 ? _v2 : _v1 ) = v; } + }; + typedef NCollection_DataMap< TGeomID, ConcaveFace > TConcaveVertex2Face; + // -------------------------------------------------------------------------- /*! * \brief Container of IDs of SOLID sub-shapes */ class Solid // sole SOLID contains all sub-shapes { - TGeomID _id; // SOLID id - bool _hasInternalFaces; + TGeomID _id; // SOLID id + bool _hasInternalFaces; + TConcaveVertex2Face _concaveVertex; // concave VERTEX -> ConcaveFace public: virtual ~Solid() {} - virtual bool Contains( TGeomID subID ) const { return true; } - virtual bool ContainsAny( const vector< TGeomID>& subIDs ) const { return true; } + virtual bool Contains( TGeomID /*subID*/ ) const { return true; } + virtual bool ContainsAny( const vector< TGeomID>& /*subIDs*/ ) const { return true; } virtual TopAbs_Orientation Orientation( const TopoDS_Shape& s ) const { return s.Orientation(); } - virtual bool IsOutsideOriented( TGeomID faceID ) const { return true; } + virtual bool IsOutsideOriented( TGeomID /*faceID*/ ) const { return true; } void SetID( TGeomID id ) { _id = id; } TGeomID ID() const { return _id; } void SetHasInternalFaces( bool has ) { _hasInternalFaces = has; } bool HasInternalFaces() const { return _hasInternalFaces; } + void SetConcave( TGeomID V, TGeomID F, TGeomID E1, TGeomID E2, TGeomID V1, TGeomID V2 ) + { _concaveVertex.Bind( V, ConcaveFace{ F, E1, E2, V1, V2 }); } + bool HasConcaveVertex() const { return !_concaveVertex.IsEmpty(); } + const ConcaveFace* GetConcave( TGeomID V ) const { return _concaveVertex.Seek( V ); } }; // -------------------------------------------------------------------------- class OneOfSolids : public Solid @@ -227,6 +254,47 @@ namespace } }; // -------------------------------------------------------------------------- + /*! + * \brief Hold a vector of TGeomID and clear it at destruction + */ + class GeomIDVecHelder + { + typedef std::vector< TGeomID > TVector; + const TVector& myVec; + bool myOwn; + + public: + GeomIDVecHelder( const TVector& idVec, bool isOwner ): myVec( idVec ), myOwn( isOwner ) {} + GeomIDVecHelder( const GeomIDVecHelder& holder ): myVec( holder.myVec ), myOwn( holder.myOwn ) + { + const_cast< bool& >( holder.myOwn ) = false; + } + ~GeomIDVecHelder() { if ( myOwn ) const_cast( myVec ).clear(); } + size_t size() const { return myVec.size(); } + TGeomID operator[]( size_t i ) const { return i < size() ? myVec[i] : theUndefID; } + bool operator==( const GeomIDVecHelder& other ) const { return myVec == other.myVec; } + bool contain( const TGeomID& id ) const { + return std::find( myVec.begin(), myVec.end(), id ) != myVec.end(); + } + TGeomID otherThan( const TGeomID& id ) const { + for ( const TGeomID& id2 : myVec ) + if ( id != id2 ) + return id2; + return theUndefID; + } + TGeomID oneCommon( const GeomIDVecHelder& other ) const { + TGeomID common = theUndefID; + for ( const TGeomID& id : myVec ) + if ( other.contain( id )) + { + if ( common != theUndefID ) + return theUndefID; + common = id; + } + return common; + } + }; + // -------------------------------------------------------------------------- /*! * \brief Geom data */ @@ -240,10 +308,13 @@ namespace TColStd_MapOfInteger _strangeEdges; // EDGEs shared by strange FACEs TGeomID _extIntFaceID; // pseudo FACE - extension of INTERNAL FACE + TopTools_DataMapOfShapeInteger _shape2NbNodes; // nb of pre-existing nodes on shapes + Controls::ElementsOnShape _edgeClassifier; Controls::ElementsOnShape _vertexClassifier; bool IsOneSolid() const { return _solidByID.size() < 2; } + GeomIDVecHelder GetSolidIDsByShapeID( const vector< TGeomID >& shapeIDs ) const; }; // -------------------------------------------------------------------------- /*! @@ -255,9 +326,10 @@ namespace mutable vector< TGeomID > _faceIDs; B_IntersectPoint(): _node(NULL) {} - void Add( const vector< TGeomID >& fIDs, const SMDS_MeshNode* n=0 ) const; - int HasCommonFace( const B_IntersectPoint * other, int avoidFace=-1 ) const; - bool IsOnFace( int faceID ) const; + bool Add( const vector< TGeomID >& fIDs, const SMDS_MeshNode* n=0 ) const; + TGeomID HasCommonFace( const B_IntersectPoint * other, TGeomID avoidFace=-1 ) const; + size_t GetCommonFaces( const B_IntersectPoint * other, TGeomID * commonFaces ) const; + bool IsOnFace( TGeomID faceID ) const; virtual ~B_IntersectPoint() {} }; // -------------------------------------------------------------------------- @@ -362,6 +434,9 @@ namespace gp_XYZ _origin; gp_Mat _invB; // inverted basis of _axes + // index shift within _nodes of nodes of a cell from the 1st node + int _nodeShift[8]; + vector< const SMDS_MeshNode* > _nodes; // mesh nodes at grid nodes vector< const F_IntersectPoint* > _gridIntP; // grid node intersection with geometry ObjectPool< E_IntersectPoint > _edgeIntPool; // intersections with EDGEs @@ -425,8 +500,11 @@ namespace const vector< TGeomID > & GetSolidIDs( TGeomID subShapeID ) const; bool IsCorrectTransition( TGeomID faceID, const Solid* solid ); bool IsBoundaryFace( TGeomID face ) const { return _geometry._boundaryFaces.Contains( face ); } - void SetOnShape( const SMDS_MeshNode* n, const F_IntersectPoint& ip, bool unset=false ); + void SetOnShape( const SMDS_MeshNode* n, const F_IntersectPoint& ip, + TopoDS_Vertex* vertex = nullptr, bool unset = false ); + void UpdateFacesOfVertex( const B_IntersectPoint& ip, const TopoDS_Vertex& vertex ); bool IsToCheckNodePos() const { return !_toAddEdges && _toCreateFaces; } + bool IsToRemoveExcessEntities() const { return !_toAddEdges; } void SetCoordinates(const vector& xCoords, const vector& yCoords, @@ -442,12 +520,14 @@ namespace */ struct CellsAroundLink { + int _iDir; int _dInd[4][3]; size_t _nbCells[3]; int _i,_j,_k; Grid* _grid; CellsAroundLink( Grid* grid, int iDir ): + _iDir( iDir ), _dInd{ {0,0,0}, {0,0,0}, {0,0,0}, {0,0,0} }, _nbCells{ grid->_coords[0].size() - 1, grid->_coords[1].size() - 1, @@ -466,7 +546,7 @@ namespace _j = j - _dInd[iL][1]; _k = k - _dInd[iL][2]; } - bool GetCell( int iL, int& i, int& j, int& k, int& cellIndex ) + bool GetCell( int iL, int& i, int& j, int& k, int& cellIndex, int& linkIndex ) { i = _i + _dInd[iL][0]; j = _j + _dInd[iL][1]; @@ -476,6 +556,7 @@ namespace k < 0 || k >= (int)_nbCells[2] ) return false; cellIndex = _grid->CellIndex( i,j,k ); + linkIndex = iL + _iDir * 4; return true; } }; @@ -603,6 +684,10 @@ namespace { return _intPoint ? _intPoint->IsOnFace( faceID ) : false; } + size_t GetCommonFaces( const B_IntersectPoint * other, TGeomID* common ) const + { + return _intPoint && other ? _intPoint->GetCommonFaces( other, common ) : 0; + } gp_Pnt Point() const { if ( const SMDS_MeshNode* n = Node() ) @@ -657,7 +742,7 @@ namespace bool _reverse; _OrientedLink( _Link* link=0, bool reverse=false ): _link(link), _reverse(reverse) {} void Reverse() { _reverse = !_reverse; } - int NbResultLinks() const { return _link->_splits.size(); } + size_t NbResultLinks() const { return _link->_splits.size(); } _OrientedLink ResultLink(int i) const { return _OrientedLink(&_link->_splits[_reverse ? NbResultLinks()-i-1 : i],_reverse); @@ -724,11 +809,11 @@ namespace { struct _Split // data of a link split { - int _linkID; // hex link ID + int _linkID; // hex link ID _Node* _nodes[2]; int _iCheckIteration; // iteration where split is tried as Hexahedron split _Link* _checkedSplit; // split set to hex links - bool _isUsed; // used in a volume + bool _isUsed; // used in a volume _Split( _Link & split, int iLink ): _linkID( iLink ), _nodes{ split._nodes[0], split._nodes[1] }, @@ -759,9 +844,13 @@ namespace // -------------------------------------------------------------------------------- struct _Face { + SMESH_Block::TShapeID _name; vector< _OrientedLink > _links; // links on GridLine's vector< _Link > _polyLinks; // links added to close a polygonal face vector< _Node* > _eIntNodes; // nodes at intersection with EDGEs + + _Face():_name( SMESH_Block::ID_NONE ) + {} bool IsPolyLink( const _OrientedLink& ol ) { return _polyLinks.empty() ? false : @@ -789,41 +878,80 @@ namespace // -------------------------------------------------------------------------------- struct _volumeDef // holder of nodes of a volume mesh element { + typedef void* _ptr; + struct _nodeDef { const SMDS_MeshNode* _node; // mesh node at hexahedron corner const B_IntersectPoint* _intPoint; + _nodeDef(): _node(0), _intPoint(0) {} _nodeDef( _Node* n ): _node( n->_node), _intPoint( n->_intPoint ) {} const SMDS_MeshNode* Node() const { return ( _intPoint && _intPoint->_node ) ? _intPoint->_node : _node; } const E_IntersectPoint* EdgeIntPnt() const { return static_cast< const E_IntersectPoint* >( _intPoint ); } + _ptr Ptr() const { return Node() ? (_ptr) Node() : (_ptr) EdgeIntPnt(); } + bool operator==(const _nodeDef& other ) const { return Ptr() == other.Ptr(); } }; + vector< _nodeDef > _nodes; vector< int > _quantities; _volumeDef* _next; // to store several _volumeDefs in a chain TGeomID _solidID; + double _size; const SMDS_MeshElement* _volume; // new volume - _volumeDef(): _next(0), _solidID(0), _volume(0) {} + vector< SMESH_Block::TShapeID > _names; // name of side a polygon originates from + + _volumeDef(): _next(0), _solidID(0), _size(0), _volume(0) {} ~_volumeDef() { delete _next; } _volumeDef( _volumeDef& other ): - _next(0), _solidID( other._solidID ), _volume( other._volume ) - { _nodes.swap( other._nodes ); _quantities.swap( other._quantities ); other._volume = 0; } + _next(0), _solidID( other._solidID ), _size( other._size ), _volume( other._volume ) + { _nodes.swap( other._nodes ); _quantities.swap( other._quantities ); other._volume = 0; + _names.swap( other._names ); } - void Set( const vector< _Node* >& nodes, const vector< int >& quant = vector< int >() ) - { _nodes.assign( nodes.begin(), nodes.end() ); _quantities = quant; } + size_t size() const { return 1 + ( _next ? _next->size() : 0 ); } // nb _volumeDef in a chain + _volumeDef* at(int index) + { return index == 0 ? this : ( _next ? _next->at(index-1) : _next ); } void Set( _Node** nodes, int nb ) { _nodes.assign( nodes, nodes + nb ); } void SetNext( _volumeDef* vd ) { if ( _next ) { _next->SetNext( vd ); } else { _next = vd; }} + + bool IsEmpty() const { return (( _nodes.empty() ) && + ( !_next || _next->IsEmpty() )); } + bool IsPolyhedron() const { return ( !_quantities.empty() || + ( _next && !_next->_quantities.empty() )); } + + + struct _linkDef: public std::pair<_ptr,_ptr> // to join polygons in removeExcessSideDivision() + { + _nodeDef _node1;//, _node2; + mutable /*const */_linkDef *_prev, *_next; + size_t _loopIndex; + + _linkDef():_prev(0), _next(0) {} + + void init( const _nodeDef& n1, const _nodeDef& n2, size_t iLoop ) + { + _node1 = n1; //_node2 = n2; + _loopIndex = iLoop; + first = n1.Ptr(); + second = n2.Ptr(); + if ( first > second ) std::swap( first, second ); + } + void setNext( _linkDef* next ) + { + _next = next; + next->_prev = this; + } + }; }; // topology of a hexahedron - int _nodeShift[8]; _Node _hexNodes [8]; _Link _hexLinks [12]; _Face _hexQuads [6]; @@ -837,7 +965,7 @@ namespace // additional nodes created at intersection points vector< _Node > _intNodes; - // nodes inside the hexahedron (at VERTEXes) + // nodes inside the hexahedron (at VERTEXes) refer to _intNodes vector< _Node* > _vIntNodes; // computed volume elements @@ -858,7 +986,7 @@ namespace Hexahedron(Grid* grid); int MakeElements(SMESH_MesherHelper& helper, const map< TGeomID, vector< TGeomID > >& edge2faceIDsMap); - void ComputeElements( const Solid* solid = 0, int solidIndex = -1 ); + void computeElements( const Solid* solid = 0, int solidIndex = -1 ); private: Hexahedron(const Hexahedron& other, size_t i, size_t j, size_t k, int cellID ); @@ -866,7 +994,7 @@ namespace void init( size_t i ); void setIJK( size_t i ); bool compute( const Solid* solid, const IsInternalFlag intFlag ); - vector< TGeomID > getSolids(); + size_t getSolids( TGeomID ids[] ); bool isCutByInternalFace( IsInternalFlag & maxFlag ); void addEdges(SMESH_MesherHelper& helper, vector< Hexahedron* >& intersectedHex, @@ -878,11 +1006,13 @@ namespace bool addIntersection( const E_IntersectPoint* ip, vector< Hexahedron* >& hexes, int ijk[], int dIJK[] ); + bool isQuadOnFace( const size_t iQuad ); bool findChain( _Node* n1, _Node* n2, _Face& quad, vector<_Node*>& chainNodes ); bool closePolygon( _Face* polygon, vector<_Node*>& chainNodes ) const; bool findChainOnEdge( const vector< _OrientedLink >& splits, const _OrientedLink& prevSplit, const _OrientedLink& avoidSplit, + const std::set< TGeomID > & concaveFaces, size_t & iS, _Face& quad, vector<_Node*>& chn); @@ -893,6 +1023,9 @@ namespace const map< TGeomID, vector< TGeomID > >& edge2faceIDsMap ); void getVolumes( vector< const SMDS_MeshElement* > & volumes ); void getBoundaryElems( vector< const SMDS_MeshElement* > & boundaryVolumes ); + void removeExcessSideDivision(const vector< Hexahedron* >& allHexa); + void removeExcessNodes(vector< Hexahedron* >& allHexa); + void preventVolumesOverlapping(); TGeomID getAnyFace() const; void cutByExtendedInternal( std::vector< Hexahedron* >& hexes, const TColStd_MapOfInteger& intEdgeIDs ); @@ -901,7 +1034,7 @@ namespace void sortVertexNodes(vector<_Node*>& nodes, _Node* curNode, TGeomID face); bool isInHole() const; bool hasStrangeEdge() const; - bool checkPolyhedronSize( bool isCutByInternalFace ) const; + bool checkPolyhedronSize( bool isCutByInternalFace, double & volSize ) const; bool addHexa (); bool addTetra(); bool addPenta(); @@ -917,6 +1050,9 @@ namespace return nodes[i]; return 0; } + bool isCorner( const _Node* node ) const { return ( node >= &_hexNodes[0] && + node - &_hexNodes[0] < 8 ); } + bool hasEdgesAround( const ConcaveFace* cf ) const; bool isImplementEdges() const { return _grid->_edgeIntPool.nbElements(); } bool isOutParam(const double uvw[3]) const; @@ -927,7 +1063,7 @@ namespace TID2Nb::iterator id2nb = id2nbMap.insert( s0 ).first; id2nb->second++; } - }; + }; // class Hexahedron #ifdef WITH_TBB // -------------------------------------------------------------------------- @@ -942,7 +1078,7 @@ namespace { for ( size_t i = r.begin(); i != r.end(); ++i ) if ( Hexahedron* hex = _hexVec[ i ] ) - hex->ComputeElements(); + hex->computeElements(); } }; // -------------------------------------------------------------------------- @@ -987,6 +1123,32 @@ namespace di = 0; } //============================================================================= + /* + * Return a vector of SOLIDS sharing given shapes + */ + GeomIDVecHelder Geometry::GetSolidIDsByShapeID( const vector< TGeomID >& theShapeIDs ) const + { + if ( theShapeIDs.size() == 1 ) + return GeomIDVecHelder( _solidIDsByShapeID[ theShapeIDs[ 0 ]], /*owner=*/false ); + + // look for an empty slot in _solidIDsByShapeID + vector< TGeomID > * resultIDs = 0; + for ( const vector< TGeomID >& vec : _solidIDsByShapeID ) + if ( vec.empty() ) + { + resultIDs = const_cast< vector< TGeomID > * >( & vec ); + break; + } + // fill in resultIDs + for ( const TGeomID& id : theShapeIDs ) + for ( const TGeomID& solid : _solidIDsByShapeID[ id ]) + { + if ( std::find( resultIDs->begin(), resultIDs->end(), solid ) == resultIDs->end() ) + resultIDs->push_back( solid ); + } + return GeomIDVecHelder( *resultIDs, /*owner=*/true ); + } + //============================================================================= /* * Remove coincident intersection points */ @@ -1036,7 +1198,7 @@ namespace switch ( ip->_transition ) { case Trans_IN: isOut = true; break; case Trans_OUT: isOut = false; break; - case Trans_TANGENT: isOut = ( prevID == 0 ); break; + case Trans_TANGENT: isOut = ( prevID != 0 ); break; case Trans_APEX: { // singularity point (apex of a cone) @@ -1059,42 +1221,47 @@ namespace return isOut ? 0 : geom._soleSolid.ID(); } - const vector< TGeomID >& solids = geom._solidIDsByShapeID[ ip->_faceIDs[ 0 ]]; + GeomIDVecHelder solids = geom.GetSolidIDsByShapeID( ip->_faceIDs ); --ip; if ( ip->_transition == Trans_INTERNAL ) return prevID; - const vector< TGeomID >& solidsBef = geom._solidIDsByShapeID[ ip->_faceIDs[ 0 ]]; + GeomIDVecHelder solidsBef = geom.GetSolidIDsByShapeID( ip->_faceIDs ); if ( ip->_transition == Trans_IN || ip->_transition == Trans_OUT ) { if ( solidsBef.size() == 1 ) - return ( solidsBef[0] == prevID ) ? 0 : solidsBef[0]; + { + if ( solidsBef[0] == prevID ) + return ip->_transition == Trans_OUT ? 0 : solidsBef[0]; + else + return solidsBef[0]; + } - return solidsBef[ solidsBef[0] == prevID ]; + if ( solids.size() == 2 ) + { + if ( solids == solidsBef ) + return theUndefID; //solids.contain( prevID ) ? solids.otherThan( prevID ) : theUndefID; + } + return solids.oneCommon( solidsBef ); } if ( solidsBef.size() == 1 ) return solidsBef[0]; - for ( size_t i = 0; i < solids.size(); ++i ) - { - vector< TGeomID >::const_iterator it = - std::find( solidsBef.begin(), solidsBef.end(), solids[i] ); - if ( it != solidsBef.end() ) - return solids[i]; - } - return 0; + return solids.oneCommon( solidsBef ); } //================================================================================ /* * Adds face IDs */ - void B_IntersectPoint::Add( const vector< TGeomID >& fIDs, + bool B_IntersectPoint::Add( const vector< TGeomID >& fIDs, const SMDS_MeshNode* n) const { + size_t prevNbF = _faceIDs.size(); + if ( _faceIDs.empty() ) _faceIDs = fIDs; else @@ -1107,12 +1274,14 @@ namespace } if ( !_node ) _node = n; + + return prevNbF < _faceIDs.size(); } //================================================================================ /* - * Returns index of a common face if any, else zero + * Return ID of a common face if any, else zero */ - int B_IntersectPoint::HasCommonFace( const B_IntersectPoint * other, int avoidFace ) const + TGeomID B_IntersectPoint::HasCommonFace( const B_IntersectPoint * other, TGeomID avoidFace ) const { if ( other ) for ( size_t i = 0; i < other->_faceIDs.size(); ++i ) @@ -1123,9 +1292,25 @@ namespace } //================================================================================ /* - * Returns \c true if \a faceID in in this->_faceIDs + * Return faces common with other point */ - bool B_IntersectPoint::IsOnFace( int faceID ) const // returns true if faceID is found + size_t B_IntersectPoint::GetCommonFaces( const B_IntersectPoint * other, TGeomID* common ) const + { + size_t nbComm = 0; + if ( !other ) + return nbComm; + if ( _faceIDs.size() > other->_faceIDs.size() ) + return other->GetCommonFaces( this, common ); + for ( const TGeomID& face : _faceIDs ) + if ( other->IsOnFace( face )) + common[ nbComm++ ] = face; + return nbComm; + } + //================================================================================ + /* + * Return \c true if \a faceID in in this->_faceIDs + */ + bool B_IntersectPoint::IsOnFace( TGeomID faceID ) const // returns true if faceID is found { vector< TGeomID >::const_iterator it = std::find( _faceIDs.begin(), _faceIDs.end(), faceID ); @@ -1340,8 +1525,7 @@ namespace } TopTools_IndexedMapOfShape faces; - if ( _toCreateFaces || isSeveralSolids ) - TopExp::MapShapes( theShapeToMesh, TopAbs_FACE, faces ); + TopExp::MapShapes( theShapeToMesh, TopAbs_FACE, faces ); // find boundary FACEs on boundary of mesh->ShapeToMesh() if ( _toCreateFaces ) @@ -1364,6 +1548,62 @@ namespace SetSolidFather( _helper->IthVertex( 1, edge ), theShapeToMesh ); } } + + // fill in _geometry._shape2NbNodes == find already meshed sub-shapes + _geometry._shape2NbNodes.Clear(); + if ( mesh->NbNodes() > 0 ) + { + for ( TopAbs_ShapeEnum type : { TopAbs_FACE, TopAbs_EDGE, TopAbs_VERTEX }) + for ( TopExp_Explorer exp( theShapeToMesh, type ); exp.More(); exp.Next() ) + { + if ( _geometry._shape2NbNodes.IsBound( exp.Current() )) + continue; + if ( SMESHDS_SubMesh* sm = mesh->GetMeshDS()->MeshElements( exp.Current() )) + if ( sm->NbNodes() > 0 ) + _geometry._shape2NbNodes.Bind( exp.Current(), sm->NbNodes() ); + } + } + + // fill in Solid::_concaveVertex + vector< TGeomID > soleSolidID( 1, _geometry._soleSolid.ID() ); + for ( int i = 1; i <= faces.Size(); ++i ) + { + const TopoDS_Face& F = TopoDS::Face( faces( i )); + TError error; + TSideVector wires = StdMeshers_FaceSide::GetFaceWires( F, *mesh, 0, error, + nullptr, nullptr, false ); + for ( StdMeshers_FaceSidePtr& wire : wires ) + { + const int nbEdges = wire->NbEdges(); + if ( nbEdges < 2 && SMESH_Algo::isDegenerated( wire->Edge(0))) + continue; + for ( int iE1 = 0; iE1 < nbEdges; ++iE1 ) + { + if ( SMESH_Algo::isDegenerated( wire->Edge( iE1 ))) continue; + int iE2 = ( iE1 + 1 ) % nbEdges; + while ( SMESH_Algo::isDegenerated( wire->Edge( iE2 ))) + iE2 = ( iE2 + 1 ) % nbEdges; + TopoDS_Vertex V = wire->FirstVertex( iE2 ); + double angle = _helper->GetAngle( wire->Edge( iE1 ), + wire->Edge( iE2 ), F, V ); + if ( angle < -5. * M_PI / 180. ) + { + TGeomID faceID = ShapeID( F ); + const vector< TGeomID > & solids = + _geometry.IsOneSolid() ? soleSolidID : GetSolidIDs( faceID ); + for ( const TGeomID & solidID : solids ) + { + Solid* solid = GetSolid( solidID ); + TGeomID V1 = ShapeID( wire->FirstVertex( iE1 )); + TGeomID V2 = ShapeID( wire->LastVertex ( iE2 )); + solid->SetConcave( ShapeID( V ), faceID, + wire->EdgeID( iE1 ), wire->EdgeID( iE2 ), V1, V2 ); + } + } + } + } + } + return; } //================================================================================ @@ -1452,8 +1692,10 @@ namespace //================================================================================ /* * Assign to geometry a node at FACE intersection + * Return a found supporting VERTEX */ - void Grid::SetOnShape( const SMDS_MeshNode* n, const F_IntersectPoint& ip, bool unset ) + void Grid::SetOnShape( const SMDS_MeshNode* n, const F_IntersectPoint& ip, + TopoDS_Vertex* vertex, bool unset ) { TopoDS_Shape s; SMESHDS_Mesh* mesh = _helper->GetMeshDS(); @@ -1465,6 +1707,8 @@ namespace { if ( unset ) mesh->UnSetNodeOnShape( n ); mesh->SetNodeOnVertex( n, TopoDS::Vertex( s )); + if ( vertex ) + *vertex = TopoDS::Vertex( s ); } else if ( _geometry._edgeClassifier.IsSatisfy( n, &s )) { @@ -1481,6 +1725,23 @@ namespace } } //================================================================================ + /* + * Fill in B_IntersectPoint::_faceIDs with all FACEs sharing a VERTEX + */ + void Grid::UpdateFacesOfVertex( const B_IntersectPoint& ip, const TopoDS_Vertex& vertex ) + { + if ( vertex.IsNull() ) + return; + std::vector< int > faceID(1); + PShapeIteratorPtr fIt = _helper->GetAncestors( vertex, *_helper->GetMesh(), + TopAbs_FACE, & _geometry._mainShape ); + while ( const TopoDS_Shape* face = fIt->next() ) + { + faceID[ 0 ] = ShapeID( *face ); + ip.Add( faceID ); + } + } + //================================================================================ /* * Initialize a classifier */ @@ -1566,8 +1827,7 @@ namespace { // state of each node of the grid relative to the geometry const size_t nbGridNodes = _coords[0].size() * _coords[1].size() * _coords[2].size(); - const TGeomID undefID = 1e+9; - vector< TGeomID > shapeIDVec( nbGridNodes, undefID ); + vector< TGeomID > shapeIDVec( nbGridNodes, theUndefID ); _nodes.resize( nbGridNodes, 0 ); _gridIntP.resize( nbGridNodes, NULL ); @@ -1629,7 +1889,9 @@ namespace gp_XYZ xyz = lineLoc + ip->_paramOnLine * lineDir; ip->_node = mesh->AddNode( xyz.X(), xyz.Y(), xyz.Z() ); ip->_indexOnLine = nodeCoord-coord0-1; - SetOnShape( ip->_node, *ip ); + TopoDS_Vertex v; + SetOnShape( ip->_node, *ip, & v ); + UpdateFacesOfVertex( *ip, v ); } // create a mesh node at ip coincident with a grid node else @@ -1666,7 +1928,7 @@ namespace { size_t nodeIndex = NodeIndex( x, y, z ); if ( !_nodes[ nodeIndex ] && - 0 < shapeIDVec[ nodeIndex ] && shapeIDVec[ nodeIndex ] < undefID ) + 0 < shapeIDVec[ nodeIndex ] && shapeIDVec[ nodeIndex ] < theUndefID ) { gp_XYZ xyz = ( _coords[0][x] * _axes[0] + _coords[1][y] * _axes[1] + @@ -1677,7 +1939,9 @@ namespace else if ( _nodes[ nodeIndex ] && _gridIntP[ nodeIndex ] /*&& !_nodes[ nodeIndex]->GetShapeID()*/ ) { - SetOnShape( _nodes[ nodeIndex ], *_gridIntP[ nodeIndex ]); + TopoDS_Vertex v; + SetOnShape( _nodes[ nodeIndex ], *_gridIntP[ nodeIndex ], & v ); + UpdateFacesOfVertex( *_gridIntP[ nodeIndex ], v ); } } @@ -2089,14 +2353,14 @@ namespace size_t i101 = i100 + dz; size_t i011 = i010 + dz; size_t i111 = i110 + dz; - _nodeShift[ SMESH_Block::ShapeIndex( SMESH_Block::ID_V000 )] = i000; - _nodeShift[ SMESH_Block::ShapeIndex( SMESH_Block::ID_V100 )] = i100; - _nodeShift[ SMESH_Block::ShapeIndex( SMESH_Block::ID_V010 )] = i010; - _nodeShift[ SMESH_Block::ShapeIndex( SMESH_Block::ID_V110 )] = i110; - _nodeShift[ SMESH_Block::ShapeIndex( SMESH_Block::ID_V001 )] = i001; - _nodeShift[ SMESH_Block::ShapeIndex( SMESH_Block::ID_V101 )] = i101; - _nodeShift[ SMESH_Block::ShapeIndex( SMESH_Block::ID_V011 )] = i011; - _nodeShift[ SMESH_Block::ShapeIndex( SMESH_Block::ID_V111 )] = i111; + grid->_nodeShift[ SMESH_Block::ShapeIndex( SMESH_Block::ID_V000 )] = i000; + grid->_nodeShift[ SMESH_Block::ShapeIndex( SMESH_Block::ID_V100 )] = i100; + grid->_nodeShift[ SMESH_Block::ShapeIndex( SMESH_Block::ID_V010 )] = i010; + grid->_nodeShift[ SMESH_Block::ShapeIndex( SMESH_Block::ID_V110 )] = i110; + grid->_nodeShift[ SMESH_Block::ShapeIndex( SMESH_Block::ID_V001 )] = i001; + grid->_nodeShift[ SMESH_Block::ShapeIndex( SMESH_Block::ID_V101 )] = i101; + grid->_nodeShift[ SMESH_Block::ShapeIndex( SMESH_Block::ID_V011 )] = i011; + grid->_nodeShift[ SMESH_Block::ShapeIndex( SMESH_Block::ID_V111 )] = i111; vector< int > idVec; // set nodes to links @@ -2112,8 +2376,10 @@ namespace int interlace[4] = { 0, 3, 1, 2 }; // to walk by links around a face: { u0, 1v, u1, 0v } for ( int faceID = SMESH_Block::ID_Fxy0; faceID <= SMESH_Block::ID_F1yz; ++faceID ) { - SMESH_Block::GetFaceEdgesIDs( faceID, idVec ); _Face& quad = _hexQuads[ SMESH_Block::ShapeIndex( faceID )]; + quad._name = (SMESH_Block::TShapeID) faceID; + + SMESH_Block::GetFaceEdgesIDs( faceID, idVec ); bool revFace = ( faceID == SMESH_Block::ID_Fxy0 || faceID == SMESH_Block::ID_Fx1z || faceID == SMESH_Block::ID_F0yz ); @@ -2141,9 +2407,6 @@ namespace _polygons.reserve(100); // to avoid reallocation; // copy topology - for ( int i = 0; i < 8; ++i ) - _nodeShift[i] = other._nodeShift[i]; - for ( int i = 0; i < 12; ++i ) { const _Link& srcLink = other._hexLinks[ i ]; @@ -2156,6 +2419,7 @@ namespace { const _Face& srcQuad = other._hexQuads[ i ]; _Face& tgtQuad = this->_hexQuads[ i ]; + tgtQuad._name = srcQuad._name; tgtQuad._links.resize(4); for ( int j = 0; j < 4; ++j ) { @@ -2167,6 +2431,8 @@ namespace } #ifdef _DEBUG_ _cellID = cellID; +#else + (void)cellID; // unused in release mode #endif } @@ -2174,8 +2440,13 @@ namespace /*! * \brief Return IDs of SOLIDs interfering with this Hexahedron */ - vector< TGeomID > Hexahedron::getSolids() + size_t Hexahedron::getSolids( TGeomID ids[] ) { + if ( _grid->_geometry.IsOneSolid() ) + { + ids[0] = _grid->GetSolid()->ID(); + return 1; + } // count intersection points belonging to each SOLID TID2Nb id2NbPoints; id2NbPoints.reserve( 3 ); @@ -2183,8 +2454,8 @@ namespace _origNodeInd = _grid->NodeIndex( _i,_j,_k ); for ( int iN = 0; iN < 8; ++iN ) { - _hexNodes[iN]._node = _grid->_nodes [ _origNodeInd + _nodeShift[iN] ]; - _hexNodes[iN]._intPoint = _grid->_gridIntP[ _origNodeInd + _nodeShift[iN] ]; + _hexNodes[iN]._node = _grid->_nodes [ _origNodeInd + _grid->_nodeShift[iN] ]; + _hexNodes[iN]._intPoint = _grid->_gridIntP[ _origNodeInd + _grid->_nodeShift[iN] ]; if ( _hexNodes[iN]._intPoint ) // intersection with a FACE { @@ -2224,12 +2495,12 @@ namespace insertAndIncrement( solidIDs[i], id2NbPoints ); } - vector< TGeomID > solids; solids.reserve( id2NbPoints.size() ); + size_t nbSolids = 0; for ( TID2Nb::iterator id2nb = id2NbPoints.begin(); id2nb != id2NbPoints.end(); ++id2nb ) if ( id2nb->second >= 3 ) - solids.push_back( id2nb->first ); + ids[ nbSolids++ ] = id2nb->first; - return solids; + return nbSolids; } //================================================================================ @@ -2332,6 +2603,7 @@ namespace { _i = i; _j = j; _k = k; + bool isCompute = solid; if ( !solid ) solid = _grid->GetSolid(); @@ -2343,8 +2615,8 @@ namespace { _hexNodes[iN]._isInternalFlags = 0; - _hexNodes[iN]._node = _grid->_nodes [ _origNodeInd + _nodeShift[iN] ]; - _hexNodes[iN]._intPoint = _grid->_gridIntP[ _origNodeInd + _nodeShift[iN] ]; + _hexNodes[iN]._node = _grid->_nodes [ _origNodeInd + _grid->_nodeShift[iN] ]; + _hexNodes[iN]._intPoint = _grid->_gridIntP[ _origNodeInd + _grid->_nodeShift[iN] ]; if ( _hexNodes[iN]._node && !solid->Contains( _hexNodes[iN]._node->GetShapeID() )) _hexNodes[iN]._node = 0; @@ -2361,6 +2633,9 @@ namespace _intNodes.clear(); _vIntNodes.clear(); + if ( !isCompute ) + return; + if ( _nbFaceIntNodes + _eIntPoints.size() > 0 && _nbFaceIntNodes + _eIntPoints.size() + _nbCornerNodes > 3) { @@ -2500,7 +2775,7 @@ namespace case 3: // at a corner { _Node& node = _hexNodes[ subEntity - SMESH_Block::ID_FirstV ]; - if ( node.Node() > 0 ) + if ( node.Node() ) { if ( node._intPoint ) node._intPoint->Add( _eIntPoints[ iP ]->_faceIDs, _eIntPoints[ iP ]->_node ); @@ -2540,7 +2815,8 @@ namespace } // loop on _eIntPoints } - else if ( 3 < _nbCornerNodes && _nbCornerNodes < 8 ) // _nbFaceIntNodes == 0 + else if (( 3 < _nbCornerNodes && _nbCornerNodes < 8 ) || // _nbFaceIntNodes == 0 + ( !_grid->_geometry.IsOneSolid() )) { _Link split; // create sub-links (_splits) of whole links @@ -2564,21 +2840,22 @@ namespace /*! * \brief Compute mesh volumes resulted from intersection of the Hexahedron */ - void Hexahedron::ComputeElements( const Solid* solid, int solidIndex ) + void Hexahedron::computeElements( const Solid* solid, int solidIndex ) { if ( !solid ) { solid = _grid->GetSolid(); if ( !_grid->_geometry.IsOneSolid() ) { - vector< TGeomID > solidIDs = getSolids(); - if ( solidIDs.size() > 1 ) + TGeomID solidIDs[20]; + size_t nbSolids = getSolids( solidIDs ); + if ( nbSolids > 1 ) { - for ( size_t i = 0; i < solidIDs.size(); ++i ) + for ( size_t i = 0; i < nbSolids; ++i ) { solid = _grid->GetSolid( solidIDs[i] ); - ComputeElements( solid, i ); - if ( !_volumeDefs._nodes.empty() && i < solidIDs.size() - 1 ) + computeElements( solid, i ); + if ( !_volumeDefs._nodes.empty() && i < nbSolids - 1 ) _volumeDefs.SetNext( new _volumeDef( _volumeDefs )); } return; @@ -2626,6 +2903,29 @@ namespace for ( int iN = 0; iN < 8; ++iN ) _hexNodes[iN]._usedInFace = 0; + if ( intFlag & IS_CUT_BY_INTERNAL_FACE && !_grid->_toAddEdges ) // Issue #19913 + preventVolumesOverlapping(); + + std::set< TGeomID > concaveFaces; // to avoid connecting nodes laying on them + + if ( solid->HasConcaveVertex() ) + { + for ( const E_IntersectPoint* ip : _eIntPoints ) + { + if ( const ConcaveFace* cf = solid->GetConcave( ip->_shapeID )) + if ( this->hasEdgesAround( cf )) + concaveFaces.insert( cf->_concaveFace ); + } + if ( concaveFaces.empty() || concaveFaces.size() * 3 < _eIntPoints.size() ) + for ( const _Node& hexNode: _hexNodes ) + { + if ( hexNode._node && hexNode._intPoint && hexNode._intPoint->_faceIDs.size() >= 3 ) + if ( const ConcaveFace* cf = solid->GetConcave( hexNode._node->GetShapeID() )) + if ( this->hasEdgesAround( cf )) + concaveFaces.insert( cf->_concaveFace ); + } + } + // Create polygons from quadrangles // -------------------------------- @@ -2633,7 +2933,8 @@ namespace vector<_Node*> chainNodes; _Face* coplanarPolyg; - bool hasEdgeIntersections = !_eIntPoints.empty(); + const bool hasEdgeIntersections = !_eIntPoints.empty(); + const bool toCheckSideDivision = isImplementEdges() || intFlag & IS_CUT_BY_INTERNAL_FACE; for ( int iF = 0; iF < 6; ++iF ) // loop on 6 sides of a hexahedron { @@ -2642,12 +2943,20 @@ namespace _polygons.resize( _polygons.size() + 1 ); _Face* polygon = &_polygons.back(); polygon->_polyLinks.reserve( 20 ); + polygon->_name = quad._name; splits.clear(); for ( int iE = 0; iE < 4; ++iE ) // loop on 4 sides of a quadrangle - for ( int iS = 0; iS < quad._links[ iE ].NbResultLinks(); ++iS ) + for ( size_t iS = 0; iS < quad._links[ iE ].NbResultLinks(); ++iS ) splits.push_back( quad._links[ iE ].ResultLink( iS )); + if ( splits.size() == 4 && + isQuadOnFace( iF )) // check if a quad on FACE is not split + { + polygon->_links.swap( splits ); + continue; // goto the next quad + } + // add splits of links to a polygon and add _polyLinks to make // polygon's boundary closed @@ -2676,6 +2985,7 @@ namespace _polygons.resize( _polygons.size() + 1 ); polygon = &_polygons.back(); polygon->_polyLinks.reserve( 20 ); + polygon->_name = quad._name; } polygon->_links.push_back( splits[ iS ] ); splits[ iS++ ]._link = 0; @@ -2691,11 +3001,12 @@ namespace n1 = split.FirstNode(); if ( n1 == n2 && n1->_intPoint && - (( n1->_intPoint->_faceIDs.size() > 1 && isImplementEdges() ) || + (( n1->_intPoint->_faceIDs.size() > 1 && toCheckSideDivision ) || ( n1->_isInternalFlags ))) { // n1 is at intersection with EDGE - if ( findChainOnEdge( splits, polygon->_links.back(), split, iS, quad, chainNodes )) + if ( findChainOnEdge( splits, polygon->_links.back(), split, concaveFaces, + iS, quad, chainNodes )) { for ( size_t i = 1; i < chainNodes.size(); ++i ) polygon->AddPolyLink( chainNodes[i-1], chainNodes[i], prevPolyg ); @@ -2843,11 +3154,12 @@ namespace } } - set usedFaceIDs; - vector< TGeomID > faces; + std::set usedFaceIDs; + std::vector< TGeomID > faces; TGeomID curFace = 0; const size_t nbQuadPolygons = _polygons.size(); E_IntersectPoint ipTmp; + std::map< TGeomID, std::vector< const B_IntersectPoint* > > tmpAddedFace; // face added to _intPoint // create polygons by making closed chains of free links size_t iPolygon = _polygons.size(); @@ -3011,7 +3323,10 @@ namespace if ( polygon._links.size() < 2 || polygon._links[0].LastNode() != polygon._links.back().FirstNode() ) - return false; // closed polygon not found -> invalid polyhedron + { + _polygons.clear(); + break; // closed polygon not found -> invalid polyhedron + } if ( polygon._links.size() == 2 ) { @@ -3106,8 +3421,11 @@ namespace for ( int iN = 0; iN < 2; ++iN ) { _Node* n = freeLinks[ iL3 ]->_link->_nodes[ iN ]; - if ( n->_intPoint ) n->_intPoint->Add( ipTmp._faceIDs ); - else n->_intPoint = &ipTmp; + bool added = false; + if ( n->_intPoint ) added = n->_intPoint->Add( ipTmp._faceIDs ); + else n->_intPoint = &ipTmp; + if ( added ) + tmpAddedFace[ ipTmp._faceIDs[0] ].push_back( n->_intPoint ); } break; } @@ -3132,8 +3450,23 @@ namespace } // end of case ( polygon._links.size() > 2 ) } // while ( nbFreeLinks > 0 ) + for ( auto & face_ip : tmpAddedFace ) + { + curFace = face_ip.first; + for ( const B_IntersectPoint* ip : face_ip.second ) + { + auto it = std::find( ip->_faceIDs.begin(), ip->_faceIDs.end(), curFace ); + if ( it != ip->_faceIDs.end() ) + ip->_faceIDs.erase( it ); + } + } + + if ( _polygons.size() < 3 ) + return false; + // check volume size - _hasTooSmall = ! checkPolyhedronSize( intFlag & IS_CUT_BY_INTERNAL_FACE ); + double volSize = 0; + _hasTooSmall = ! checkPolyhedronSize( intFlag & IS_CUT_BY_INTERNAL_FACE, volSize ); for ( size_t i = 0; i < 8; ++i ) if ( _hexNodes[ i ]._intPoint == &ipTmp ) @@ -3142,11 +3475,45 @@ namespace if ( _hasTooSmall ) return false; // too small volume + + // Try to find out names of no-name polygons (issue # 19887) + if ( _grid->IsToRemoveExcessEntities() && _polygons.back()._name == SMESH_Block::ID_NONE ) + { + gp_XYZ uvwCenter = + 0.5 * ( _grid->_coords[0][_i] + _grid->_coords[0][_i+1] ) * _grid->_axes[0] + + 0.5 * ( _grid->_coords[1][_j] + _grid->_coords[1][_j+1] ) * _grid->_axes[1] + + 0.5 * ( _grid->_coords[2][_k] + _grid->_coords[2][_k+1] ) * _grid->_axes[2]; + for ( size_t i = _polygons.size() - 1; _polygons[i]._name == SMESH_Block::ID_NONE; --i ) + { + _Face& face = _polygons[ i ]; + Bnd_Box bb; + gp_Pnt uvw; + for ( size_t iL = 0; iL < face._links.size(); ++iL ) + { + _Node* n = face._links[ iL ].FirstNode(); + gp_XYZ p = SMESH_NodeXYZ( n->Node() ); + _grid->ComputeUVW( p, uvw.ChangeCoord().ChangeData() ); + bb.Add( uvw ); + } + gp_Pnt pMin = bb.CornerMin(); + if ( bb.IsXThin( _grid->_tol )) + face._name = pMin.X() < uvwCenter.X() ? SMESH_Block::ID_F0yz : SMESH_Block::ID_F1yz; + else if ( bb.IsYThin( _grid->_tol )) + face._name = pMin.Y() < uvwCenter.Y() ? SMESH_Block::ID_Fx0z : SMESH_Block::ID_Fx1z; + else if ( bb.IsZThin( _grid->_tol )) + face._name = pMin.Z() < uvwCenter.Z() ? SMESH_Block::ID_Fxy0 : SMESH_Block::ID_Fxy1; + } + } + + _volumeDefs._nodes.clear(); + _volumeDefs._quantities.clear(); + _volumeDefs._names.clear(); + // create a classic cell if possible int nbPolygons = 0; for ( size_t iF = 0; iF < _polygons.size(); ++iF ) - nbPolygons += (_polygons[ iF ]._links.size() > 0 ); + nbPolygons += (_polygons[ iF ]._links.size() > 2 ); //const int nbNodes = _nbCornerNodes + nbIntersections; int nbNodes = 0; @@ -3162,19 +3529,18 @@ namespace else if ( nbNodes == 5 && nbPolygons == 5 ) isClassicElem = addPyra (); if ( !isClassicElem ) { - _volumeDefs._nodes.clear(); - _volumeDefs._quantities.clear(); - for ( size_t iF = 0; iF < _polygons.size(); ++iF ) { const size_t nbLinks = _polygons[ iF ]._links.size(); - if ( nbLinks == 0 ) continue; + if ( nbLinks < 3 ) continue; _volumeDefs._quantities.push_back( nbLinks ); + _volumeDefs._names.push_back( _polygons[ iF ]._name ); for ( size_t iL = 0; iL < nbLinks; ++iL ) _volumeDefs._nodes.push_back( _polygons[ iF ]._links[ iL ].FirstNode() ); } } _volumeDefs._solidID = solid->ID(); + _volumeDefs._size = volSize; return !_volumeDefs._nodes.empty(); } @@ -3193,7 +3559,7 @@ namespace int nbIntHex = 0; // set intersection nodes from GridLine's to links of allHexa - int i,j,k, cellIndex; + int i,j,k, cellIndex, iLink; for ( int iDir = 0; iDir < 3; ++iDir ) { // loop on GridLine's parallel to iDir @@ -3210,7 +3576,7 @@ namespace fourCells.Init( lineInd.I(), lineInd.J(), lineInd.K() ); for ( int iL = 0; iL < 4; ++iL ) // loop on 4 cells sharing a link { - if ( !fourCells.GetCell( iL, i,j,k, cellIndex )) + if ( !fourCells.GetCell( iL, i,j,k, cellIndex, iLink )) continue; Hexahedron *& hex = allHexa[ cellIndex ]; if ( !hex) @@ -3218,7 +3584,6 @@ namespace hex = new Hexahedron( *this, i, j, k, cellIndex ); ++nbIntHex; } - const int iLink = iL + iDir * 4; hex->_hexLinks[iLink]._fIntPoints.push_back( &(*ip) ); hex->_nbFaceIntNodes += bool( ip->_node ); } @@ -3231,6 +3596,7 @@ namespace // add not split hexahedra to the mesh int nbAdded = 0; + TGeomID solidIDs[20]; vector< Hexahedron* > intHexa; intHexa.reserve( nbIntHex ); vector< const SMDS_MeshElement* > boundaryVolumes; boundaryVolumes.reserve( nbIntHex * 1.1 ); for ( size_t i = 0; i < allHexa.size(); ++i ) @@ -3240,7 +3606,9 @@ namespace if ( hex ) // split hexahedron { intHexa.push_back( hex ); - if ( hex->_nbFaceIntNodes > 0 || hex->_eIntPoints.size() > 0 ) + if ( hex->_nbFaceIntNodes > 0 || + hex->_eIntPoints.size() > 0 || + hex->getSolids( solidIDs ) > 1 ) continue; // treat intersected hex later in parallel this->init( hex->_i, hex->_j, hex->_k ); } @@ -3266,7 +3634,8 @@ namespace } else { - solidID = getSolids()[0]; + getSolids( solidIDs ); + solidID = solidIDs[0]; } mesh->SetMeshElementOnShape( el, solidID ); ++nbAdded; @@ -3280,28 +3649,39 @@ namespace } else if ( _nbCornerNodes > 3 && !hex ) { - // all intersection of hex with geometry are at grid nodes + // all intersections of hex with geometry are at grid nodes hex = new Hexahedron( *this, _i, _j, _k, i ); intHexa.push_back( hex ); } } - // add elements resulted from hexadron intersection + // compute definitions of volumes resulted from hexadron intersection #ifdef WITH_TBB tbb::parallel_for ( tbb::blocked_range( 0, intHexa.size() ), ParallelHexahedron( intHexa ), - tbb::simple_partitioner()); // ComputeElements() is called here + tbb::simple_partitioner()); // computeElements() is called here +#else for ( size_t i = 0; i < intHexa.size(); ++i ) if ( Hexahedron * hex = intHexa[ i ] ) - nbAdded += hex->addVolumes( helper ); -#else + hex->computeElements(); +#endif + + // simplify polyhedrons + if ( _grid->IsToRemoveExcessEntities() ) + { + for ( size_t i = 0; i < intHexa.size(); ++i ) + if ( Hexahedron * hex = intHexa[ i ] ) + hex->removeExcessSideDivision( allHexa ); + + for ( size_t i = 0; i < intHexa.size(); ++i ) + if ( Hexahedron * hex = intHexa[ i ] ) + hex->removeExcessNodes( allHexa ); + } + + // add volumes for ( size_t i = 0; i < intHexa.size(); ++i ) if ( Hexahedron * hex = intHexa[ i ] ) - { - hex->ComputeElements(); nbAdded += hex->addVolumes( helper ); - } -#endif // fill boundaryVolumes with volumes neighboring too small skipped volumes if ( _grid->_toCreateFaces ) @@ -3311,6 +3691,22 @@ namespace hex->getBoundaryElems( boundaryVolumes ); } + // merge nodes on outer sub-shapes with pre-existing ones + TopTools_DataMapIteratorOfDataMapOfShapeInteger s2nIt( _grid->_geometry._shape2NbNodes ); + for ( ; s2nIt.More(); s2nIt.Next() ) + if ( s2nIt.Value() > 0 ) + if ( SMESHDS_SubMesh* sm = mesh->MeshElements( s2nIt.Key() )) + { + TIDSortedNodeSet smNodes( SMDS_MeshElement::iterator( sm->GetNodes() ), + SMDS_MeshElement::iterator() ); + SMESH_MeshEditor::TListOfListOfNodes equalNodes; + SMESH_MeshEditor editor( helper.GetMesh() ); + editor.FindCoincidentNodes( smNodes, 10 * _grid->_tol, equalNodes, + /*SeparateCornersAndMedium =*/ false); + if ((int) equalNodes.size() <= s2nIt.Value() ) + editor.MergeNodes( equalNodes ); + } + // create boundary mesh faces addFaces( helper, boundaryVolumes ); @@ -3386,7 +3782,7 @@ namespace continue; // perform intersection - E_IntersectPoint* eip, *vip; + E_IntersectPoint* eip, *vip = 0; for ( int iDirZ = 0; iDirZ < 3; ++iDirZ ) { GridPlanes& planes = pln[ iDirZ ]; @@ -3424,6 +3820,7 @@ namespace ip._point = p1; ip._shapeID = _grid->ShapeID( v1 ); vip = _grid->Add( ip ); + _grid->UpdateFacesOfVertex( *vip, v1 ); if ( isInternal ) vip->_faceIDs.push_back( _grid->PseudoIntExtFaceID() ); if ( !addIntersection( vip, hexes, ijk, d000 )) @@ -3484,9 +3881,12 @@ namespace ijk[ iDirZ ] = iZ1; bool sameV = ( v1.IsSame( v2 )); if ( !sameV ) + { vip = _grid->Add( ip ); - if ( isInternal && !sameV ) - vip->_faceIDs.push_back( _grid->PseudoIntExtFaceID() ); + _grid->UpdateFacesOfVertex( *vip, v2 ); + if ( isInternal ) + vip->_faceIDs.push_back( _grid->PseudoIntExtFaceID() ); + } if ( !addIntersection( vip, hexes, ijk, d000 ) && !sameV ) _grid->Remove( vip ); ip._shapeID = edgeID; @@ -3596,7 +3996,8 @@ namespace int i = ! ( u < _grid->_tol ); // [0,1] int iN = link._nodes[ i ] - hex->_hexNodes; // [0-7] - const F_IntersectPoint * & ip = _grid->_gridIntP[ hex->_origNodeInd + _nodeShift[iN] ]; + const F_IntersectPoint * & ip = _grid->_gridIntP[ hex->_origNodeInd + + _grid->_nodeShift[iN] ]; if ( !ip ) { ip = _grid->_extIntPool.getNew(); @@ -3624,13 +4025,12 @@ namespace int i,j,k, cellIndex; for ( int iC = 0; iC < 4; ++iC ) // loop on 4 cells sharing the link { - if ( !fourCells.GetCell( iC, i,j,k, cellIndex )) + if ( !fourCells.GetCell( iC, i,j,k, cellIndex, iLink )) continue; Hexahedron * h = hexes[ cellIndex ]; if ( !h ) h = hexes[ cellIndex ] = new Hexahedron( *this, i, j, k, cellIndex ); - const int iL = iC + iDir * 4; - h->_hexLinks[iL]._fIntPoints.push_back( ip ); + h->_hexLinks[iLink]._fIntPoints.push_back( ip ); h->_nbFaceIntNodes++; //isCut = true; } @@ -3690,8 +4090,8 @@ namespace for ( size_t iN = 0; iN < 8; ++iN ) // check corners { - _hexNodes[iN]._node = _grid->_nodes [ _origNodeInd + _nodeShift[iN] ]; - _hexNodes[iN]._intPoint = _grid->_gridIntP[ _origNodeInd + _nodeShift[iN] ]; + _hexNodes[iN]._node = _grid->_nodes [ _origNodeInd + _grid->_nodeShift[iN] ]; + _hexNodes[iN]._intPoint = _grid->_gridIntP[ _origNodeInd + _grid->_nodeShift[iN] ]; if ( _hexNodes[iN]._intPoint ) for ( size_t iF = 0; iF < _hexNodes[iN]._intPoint->_faceIDs.size(); ++iF ) { @@ -3859,6 +4259,55 @@ namespace return added; } //================================================================================ + /*! + * \brief Check if a hexahedron facet lies on a FACE + * Also return true if the facet does not interfere with any FACE + */ + bool Hexahedron::isQuadOnFace( const size_t iQuad ) + { + _Face& quad = _hexQuads[ iQuad ] ; + + int nbGridNodesInt = 0; // nb FACE intersections at grid nodes + int nbNoGeomNodes = 0; + for ( int iE = 0; iE < 4; ++iE ) + { + nbNoGeomNodes = ( !quad._links[ iE ].FirstNode()->_intPoint && + quad._links[ iE ].NbResultLinks() == 1 ); + nbGridNodesInt += + ( quad._links[ iE ].FirstNode()->_intPoint && + quad._links[ iE ].NbResultLinks() == 1 && + quad._links[ iE ].ResultLink( 0 ).FirstNode() == quad._links[ iE ].FirstNode() && + quad._links[ iE ].ResultLink( 0 ).LastNode() == quad._links[ iE ].LastNode() ); + } + if ( nbNoGeomNodes == 4 ) + return true; + + if ( nbGridNodesInt == 4 ) // all quad nodes are at FACE intersection + { + size_t iEmin = 0, minNbFaces = 1000; + for ( int iE = 0; iE < 4; ++iE ) // look for a node with min nb FACEs + { + size_t nbFaces = quad._links[ iE ].FirstNode()->faces().size(); + if ( minNbFaces > nbFaces ) + { + iEmin = iE; + minNbFaces = nbFaces; + } + } + // check if there is a FACE passing through all 4 nodes + for ( const TGeomID& faceID : quad._links[ iEmin ].FirstNode()->faces() ) + { + bool allNodesAtFace = true; + for ( size_t iE = 0; iE < 4 && allNodesAtFace; ++iE ) + allNodesAtFace = ( iE == iEmin || + quad._links[ iE ].FirstNode()->IsOnFace( faceID )); + if ( allNodesAtFace ) // quad if on faceID + return true; + } + } + return false; + } + //================================================================================ /*! * \brief Finds nodes at a path from one node to another via intersections with EDGEs */ @@ -3909,8 +4358,8 @@ namespace return false; vector< _OrientedLink > newLinks; // find a node lying on the same FACE as the last one - _Node* node = polygon->_links.back().LastNode(); - int avoidFace = node->IsLinked( polygon->_links.back().FirstNode()->_intPoint ); + _Node* node = polygon->_links.back().LastNode(); + TGeomID avoidFace = node->IsLinked( polygon->_links.back().FirstNode()->_intPoint ); for ( i = nbLinks - 2; i >= 0; --i ) if ( node->IsLinked( polygon->_links[i].FirstNode()->_intPoint, avoidFace )) break; @@ -3953,19 +4402,20 @@ namespace bool Hexahedron::findChainOnEdge( const vector< _OrientedLink >& splits, const _OrientedLink& prevSplit, const _OrientedLink& avoidSplit, + const std::set< TGeomID > & concaveFaces, size_t & iS, _Face& quad, vector<_Node*>& chn ) { _Node* pn1 = prevSplit.FirstNode(); - _Node* pn2 = prevSplit.LastNode(); - int avoidFace = pn1->IsLinked( pn2->_intPoint ); // FACE under the quad + _Node* pn2 = prevSplit.LastNode(); // pn2 is on EDGE, if not on INTERNAL FACE + _Node* an3 = avoidSplit.LastNode(); + TGeomID avoidFace = pn1->IsLinked( pn2->_intPoint ); // FACE under the quad if ( avoidFace < 1 && pn1->_intPoint ) return false; - _Node* n = 0, *stopNode = avoidSplit.LastNode(); - chn.clear(); + if ( !quad._eIntNodes.empty() ) // connect pn2 with EDGE intersections { chn.push_back( pn2 ); @@ -3986,28 +4436,82 @@ namespace pn2 = chn.back(); } - int i; - for ( i = splits.size()-1; i >= 0; --i ) // connect new pn2 (at _eIntNodes) with a split + _Node* n = 0, *stopNode = avoidSplit.LastNode(); + + if ( pn2 == prevSplit.LastNode() && // pn2 is at avoidSplit.FirstNode() + !isCorner( stopNode )) // stopNode is in the middle of a _hexLinks + { + // move stopNode to a _hexNodes + for ( int iE = 0; iE < 4; ++iE ) // loop on 4 sides of a quadrangle + for ( size_t iL = 0; iL < quad._links[ iE ].NbResultLinks(); ++iL ) + { + const _Link* sideSplit = & quad._links[ iE ]._link->_splits[ iL ]; + if ( sideSplit == avoidSplit._link ) + { + if ( quad._links[ iE ].LastNode()->Node() ) + stopNode = quad._links[ iE ].LastNode(); + iE = 4; + break; + } + } + } + + // connect pn2 (probably new, at _eIntNodes) with a split + + int i, iConn; + size_t nbCommon; + TGeomID commonFaces[20]; + _Node* nPrev = nullptr; + for ( i = splits.size()-1; i >= 0; --i ) { if ( !splits[i] ) continue; - n = splits[i].LastNode(); - if ( n == stopNode ) - break; - if (( n != pn1 ) && - ( n->IsLinked( pn2->_intPoint, avoidFace )) && - ( !avoidFace || n->IsOnFace( avoidFace ))) - break; + bool stop = false; + for ( int is1st = 0; is1st < 2; ++is1st ) + { + _Node* nConn = is1st ? splits[i].FirstNode() : splits[i].LastNode(); + if ( nConn == nPrev ) + { + if ( n == nConn ) + iConn = i; + continue; + } + nPrev = nConn; + if (( stop = ( nConn == stopNode ))) + break; + // find a FACE connecting nConn with pn2 but not with an3 + if (( nConn != pn1 ) && + ( nConn->_intPoint && !nConn->_intPoint->_faceIDs.empty() ) && + ( nbCommon = nConn->GetCommonFaces( pn2->_intPoint, commonFaces ))) + { + bool a3Coonect = true; + for ( size_t iF = 0; iF < nbCommon && a3Coonect; ++iF ) + a3Coonect = an3->IsOnFace( commonFaces[ iF ]) || concaveFaces.count( commonFaces[ iF ]); + if ( a3Coonect ) + continue; - n = splits[i].FirstNode(); - if ( n == stopNode ) - break; - if (( n->IsLinked( pn2->_intPoint, avoidFace )) && - ( !avoidFace || n->IsOnFace( avoidFace ))) + if ( !n ) + { + n = nConn; + iConn = i + !is1st; + } + if ( nbCommon > 1 ) // nConn is linked with pn2 by an EDGE + { + n = nConn; + iConn = i + !is1st; + stop = true; + break; + } + } + } + if ( stop ) + { + i = iConn; break; - n = 0; + } } + if ( n && n != stopNode ) { if ( chn.empty() ) @@ -4019,8 +4523,8 @@ namespace else if ( !chn.empty() && chn.back()->_isInternalFlags ) { // INTERNAL FACE partially cuts the quad - for ( int i = chn.size() - 2; i >= 0; --i ) - chn.push_back( chn[ i ]); + for ( int ip = chn.size() - 2; ip >= 0; --ip ) + chn.push_back( chn[ ip ]); return true; } return false; @@ -4291,18 +4795,34 @@ namespace mesh->SetNodeOnEdge( nodes[iN], shapeID ); } else if ( toCheckNodePos && - !nodes[iN]->isMarked() && + !nodes[iN]->isMarked() && _grid->ShapeType( nodes[iN]->GetShapeID() ) == TopAbs_FACE ) { - _grid->SetOnShape( nodes[iN], noIntPnt, /*unset=*/true ); + _grid->SetOnShape( nodes[iN], noIntPnt, /*v=*/nullptr,/*unset=*/true ); nodes[iN]->setIsMarked( true ); } - } + } // loop to get nodes const SMDS_MeshElement* v = 0; if ( !volDef->_quantities.empty() ) { v = helper.AddPolyhedralVolume( nodes, volDef->_quantities ); + volDef->_size = SMDS_VolumeTool( v ).GetSize(); + if ( volDef->_size < 0 ) // invalid polyhedron + { + if ( ! SMESH_MeshEditor( helper.GetMesh() ).Reorient( v ) || // try to fix + SMDS_VolumeTool( v ).GetSize() < 0 ) + { + helper.GetMeshDS()->RemoveFreeElement( v, /*sm=*/nullptr, /*fromGroups=*/false ); + v = nullptr; + //_hasTooSmall = true; +#ifdef _DEBUG_ + std::cout << "Remove INVALID polyhedron, _cellID = " << _cellID + << " ijk = ( " << _i << " " << _j << " " << _k << " ) " + << " solid " << volDef->_solidID << std::endl; +#endif + } + } } else { @@ -4319,10 +4839,46 @@ namespace break; } } - if (( volDef->_volume = v )) + volDef->_volume = v; + nbAdded += bool( v ); + + } // loop on _volumeDefs chain + + // avoid creating overlapping volumes (bos #24052) + if ( nbAdded > 1 ) + { + double sumSize = 0, maxSize = 0; + _volumeDef* maxSizeDef = nullptr; + for ( _volumeDef* volDef = &_volumeDefs; volDef; volDef = volDef->_next ) { - helper.GetMeshDS()->SetMeshElementOnShape( v, volDef->_solidID ); - ++nbAdded; + if ( !volDef->_volume ) + continue; + sumSize += volDef->_size; + if ( volDef->_size > maxSize ) + { + maxSize = volDef->_size; + maxSizeDef = volDef; + } + } + if ( sumSize > _sideLength[0] * _sideLength[1] * _sideLength[2] * 1.05 ) + { + for ( _volumeDef* volDef = &_volumeDefs; volDef; volDef = volDef->_next ) + if ( volDef != maxSizeDef && volDef->_volume ) + { + helper.GetMeshDS()->RemoveFreeElement( volDef->_volume, /*sm=*/nullptr, + /*fromGroups=*/false ); + volDef->_volume = nullptr; + //volDef->_nodes.clear(); + --nbAdded; + } + } + } + + for ( _volumeDef* volDef = &_volumeDefs; volDef; volDef = volDef->_next ) + { + if ( volDef->_volume ) + { + helper.GetMeshDS()->SetMeshElementOnShape( volDef->_volume, volDef->_solidID ); } } @@ -4361,10 +4917,13 @@ namespace { curIntPnt._paramOnLine = coords[ ijk[ iDir ]] - coords[0] + _grid->_tol; const GridLine& line = _grid->_lines[ iDir ][ lineIndex[ iL ]]; - multiset< F_IntersectPoint >::const_iterator ip = - line._intPoints.upper_bound( curIntPnt ); - --ip; - firstIntPnt = &(*ip); + if ( !line._intPoints.empty() ) + { + multiset< F_IntersectPoint >::const_iterator ip = + line._intPoints.upper_bound( curIntPnt ); + --ip; + firstIntPnt = &(*ip); + } } else if ( !link._fIntPoints.empty() ) { @@ -4422,8 +4981,10 @@ namespace /*! * \brief Return true if a polyhedron passes _sizeThreshold criterion */ - bool Hexahedron::checkPolyhedronSize( bool cutByInternalFace ) const + bool Hexahedron::checkPolyhedronSize( bool cutByInternalFace, double & volume) const { + volume = 0; + if ( cutByInternalFace && !_grid->_toUseThresholdForInternalFaces ) { // check if any polygon fully lies on shared/internal FACEs @@ -4443,10 +5004,6 @@ namespace return true; } } - if ( this->hasStrangeEdge() ) - return true; - - double volume = 0; for ( size_t iP = 0; iP < _polygons.size(); ++iP ) { const _Face& polygon = _polygons[iP]; @@ -4464,6 +5021,9 @@ namespace } volume /= 6; + if ( this->hasStrangeEdge() && volume > 1e-13 ) + return true; + double initVolume = _sideLength[0] * _sideLength[1] * _sideLength[2]; return volume > initVolume / _grid->_sizeThreshold; @@ -4631,6 +5191,63 @@ namespace return false; } //================================================================================ + /*! + * \brief Return true if there are _eIntPoints at EDGEs forming a concave corner + */ + bool Hexahedron::hasEdgesAround( const ConcaveFace* cf ) const + { + int nbEdges = 0; + ConcaveFace foundGeomHolder; + for ( const E_IntersectPoint* ip : _eIntPoints ) + { + if ( cf->HasEdge( ip->_shapeID )) + { + if ( ++nbEdges == 2 ) + return true; + foundGeomHolder.SetEdge( ip->_shapeID ); + } + else if ( ip->_faceIDs.size() >= 3 ) + { + const TGeomID & vID = ip->_shapeID; + if ( cf->HasVertex( vID ) && !foundGeomHolder.HasVertex( vID )) + { + if ( ++nbEdges == 2 ) + return true; + foundGeomHolder.SetVertex( vID ); + } + } + } + + for ( const _Node& hexNode: _hexNodes ) + { + if ( !hexNode._node || !hexNode._intPoint ) + continue; + const B_IntersectPoint* ip = hexNode._intPoint; + if ( ip->_faceIDs.size() == 2 ) // EDGE + { + TGeomID edgeID = hexNode._node->GetShapeID(); + if ( cf->HasEdge( edgeID ) && !foundGeomHolder.HasEdge( edgeID )) + { + foundGeomHolder.SetEdge( edgeID ); + if ( ++nbEdges == 2 ) + return true; + } + } + else if ( ip->_faceIDs.size() >= 3 ) // VERTEX + { + TGeomID vID = hexNode._node->GetShapeID(); + if ( cf->HasVertex( vID ) && !foundGeomHolder.HasVertex( vID )) + { + if ( ++nbEdges == 2 ) + return true; + foundGeomHolder.SetVertex( vID ); + } + } + } + + return false; + } + //================================================================================ /*! * \brief Dump a link and return \c false */ @@ -4641,6 +5258,8 @@ namespace cout << "BUG: not shared link. IKJ = ( "<< _i << " " << _j << " " << _k << " )" << endl << "n1 (" << p1.X() << ", "<< p1.Y() << ", "<< p1.Z() << " )" << endl << "n2 (" << p2.X() << ", "<< p2.Y() << ", "<< p2.Z() << " )" << endl; +#else + (void)link; // unused in release mode #endif return false; } @@ -4658,6 +5277,89 @@ namespace ( _grid->_coords[2][ _k+1 ] + _grid->_tol < uvw[2] )); } //================================================================================ + /*! + * \brief Find existing triangulation of a polygon + */ + int findExistingTriangulation( const SMDS_MeshElement* polygon, + //const SMDS_Mesh* mesh, + std::vector< const SMDS_MeshNode* >& nodes ) + { + int nbSplits = 0; + nodes.clear(); + std::vector twoNodes(2); + std::vector foundFaces; foundFaces.reserve(10); + std::set< const SMDS_MeshElement * > avoidFaces; avoidFaces.insert( polygon ); + + const int nbPolyNodes = polygon->NbCornerNodes(); + twoNodes[1] = polygon->GetNode( nbPolyNodes - 1 ); + for ( int iN = 0; iN < nbPolyNodes; ++iN ) // loop on border links of polygon + { + twoNodes[0] = polygon->GetNode( iN ); + + int nbFaces = SMDS_Mesh::GetElementsByNodes( twoNodes, foundFaces, SMDSAbs_Face ); + int nbOkFaces = 0; + for ( int iF = 0; iF < nbFaces; ++iF ) // keep faces lying over polygon + { + if ( avoidFaces.count( foundFaces[ iF ])) + continue; + int i, nbFaceNodes = foundFaces[ iF ]->NbCornerNodes(); + for ( i = 0; i < nbFaceNodes; ++i ) + { + const SMDS_MeshNode* n = foundFaces[ iF ]->GetNode( i ); + bool isCommonNode = ( n == twoNodes[0] || + n == twoNodes[1] || + polygon->GetNodeIndex( n ) >= 0 ); + if ( !isCommonNode ) + break; + } + if ( i == nbFaceNodes ) // all nodes of foundFaces[iF] are shared with polygon + if ( nbOkFaces++ != iF ) + foundFaces[ nbOkFaces-1 ] = foundFaces[ iF ]; + } + if ( nbOkFaces > 0 ) + { + int iFaceSelected = 0; + if ( nbOkFaces > 1 ) // select a face with minimal distance from polygon + { + double minDist = Precision::Infinite(); + for ( int iF = 0; iF < nbOkFaces; ++iF ) + { + int i, nbFaceNodes = foundFaces[ iF ]->NbCornerNodes(); + gp_XYZ gc = SMESH_NodeXYZ( foundFaces[ iF ]->GetNode( 0 )); + for ( i = 1; i < nbFaceNodes; ++i ) + gc += SMESH_NodeXYZ( foundFaces[ iF ]->GetNode( i )); + gc /= nbFaceNodes; + + double dist = SMESH_MeshAlgos::GetDistance( polygon, gc ); + if ( dist < minDist ) + { + minDist = dist; + iFaceSelected = iF; + } + } + } + if ( foundFaces[ iFaceSelected ]->NbCornerNodes() != 3 ) + return 0; + nodes.insert( nodes.end(), + foundFaces[ iFaceSelected ]->begin_nodes(), + foundFaces[ iFaceSelected ]->end_nodes()); + if ( !SMESH_MeshAlgos::IsRightOrder( foundFaces[ iFaceSelected ], + twoNodes[0], twoNodes[1] )) + { + // reverse just added nodes + std::reverse( nodes.end() - 3, nodes.end() ); + } + avoidFaces.insert( foundFaces[ iFaceSelected ]); + nbSplits++; + } + + twoNodes[1] = twoNodes[0]; + + } // loop on polygon nodes + + return nbSplits; + } + //================================================================================ /*! * \brief Divide a polygon into triangles and modify accordingly an adjacent polyhedron */ @@ -4671,7 +5373,12 @@ namespace const bool reinitVolume) { SMESH_MeshAlgos::Triangulate divider(/*optimize=*/false); - int nbTrias = divider.GetTriangles( polygon, face.myNodes ); + bool triangulationExist = false; + int nbTrias = findExistingTriangulation( polygon, face.myNodes ); + if ( nbTrias > 0 ) + triangulationExist = true; + else + nbTrias = divider.GetTriangles( polygon, face.myNodes ); face.myNodes.resize( nbTrias * 3 ); SMESH_MeshEditor::ElemFeatures newVolumeDef; @@ -4692,8 +5399,11 @@ namespace face.myNodes.begin(), face.myNodes.begin() + 3 ); meshDS->RemoveFreeElement( polygon, 0, false ); - newTriangle = meshDS->AddFace( face.myNodes[0], face.myNodes[1], face.myNodes[2] ); - meshDS->SetMeshElementOnShape( newTriangle, faceID ); + if ( !triangulationExist ) + { + newTriangle = meshDS->AddFace( face.myNodes[0], face.myNodes[1], face.myNodes[2] ); + meshDS->SetMeshElementOnShape( newTriangle, faceID ); + } } else { @@ -4710,8 +5420,11 @@ namespace newVolumeDef.myNodes.insert( newVolumeDef.myNodes.end(), face.myNodes.begin() + iN, face.myNodes.begin() + iN + 3 ); - newTriangle = meshDS->AddFace( face.myNodes[iN], face.myNodes[iN+1], face.myNodes[iN+2] ); - meshDS->SetMeshElementOnShape( newTriangle, faceID ); + if ( !triangulationExist ) + { + newTriangle = meshDS->AddFace( face.myNodes[iN], face.myNodes[iN+1], face.myNodes[iN+2] ); + meshDS->SetMeshElementOnShape( newTriangle, faceID ); + } } meshDS->RemoveFreeElement( volume.Element(), 0, false ); @@ -4726,6 +5439,36 @@ namespace return; } //================================================================================ + /*! + * \brief Look for a FACE supporting all given nodes made on EDGEs and VERTEXes + */ + TGeomID findCommonFace( const std::vector< const SMDS_MeshNode* > & nn, + const SMESH_Mesh* mesh ) + { + TGeomID faceID = 0; + TGeomID shapeIDs[20]; + for ( size_t iN = 0; iN < nn.size(); ++iN ) + shapeIDs[ iN ] = nn[ iN ]->GetShapeID(); + + SMESH_subMesh* sm = mesh->GetSubMeshContaining( shapeIDs[ 0 ]); + for ( const SMESH_subMesh * smFace : sm->GetAncestors() ) + { + if ( smFace->GetSubShape().ShapeType() != TopAbs_FACE ) + continue; + + faceID = smFace->GetId(); + + for ( size_t iN = 1; iN < nn.size() && faceID; ++iN ) + { + if ( !smFace->DependsOn( shapeIDs[ iN ])) + faceID = 0; + } + if ( faceID > 0 ) + break; + } + return faceID; + } + //================================================================================ /*! * \brief Create mesh faces at free facets */ @@ -4758,12 +5501,14 @@ namespace bndFacets.clear(); for ( int iF = 0, n = vTool.NbFaces(); iF < n; iF++ ) { - bool isBoundary = vTool.IsFreeFace( iF ); + const SMDS_MeshElement* otherVol; + bool isBoundary = vTool.IsFreeFace( iF, &otherVol ); if ( isBoundary ) { bndFacets.push_back( iF ); } - else if ( hasInternal ) + else if (( hasInternal ) || + ( !_grid->IsSolid( otherVol->GetShapeID() ))) { // check if all nodes are on internal/shared FACEs isBoundary = true; @@ -4807,18 +5552,8 @@ namespace if ( nn[ iN ]->GetPosition()->GetDim() == 2 ) faceID = nn[ iN ]->GetShapeID(); } - for ( size_t iN = 0; iN < nbFaceNodes && !faceID; ++iN ) - { - // look for a father FACE of EDGEs and VERTEXes - const TopoDS_Shape& s1 = _grid->Shape( nn[ iN ]->GetShapeID() ); - const TopoDS_Shape& s2 = _grid->Shape( nn[ iN+1 ]->GetShapeID() ); - if ( s1 != s2 && s1.ShapeType() == TopAbs_EDGE && s2.ShapeType() == TopAbs_EDGE ) - { - TopoDS_Shape f = helper.GetCommonAncestor( s1, s2, *helper.GetMesh(), TopAbs_FACE ); - if ( !f.IsNull() ) - faceID = _grid->ShapeID( f ); - } - } + if ( faceID == 0 ) + faceID = findCommonFace( face.myNodes, helper.GetMesh() ); bool toCheckFace = faceID && (( !isBoundary ) || ( hasInternal && _grid->_toUseThresholdForInternalFaces )); @@ -4831,12 +5566,15 @@ namespace if ( subID != faceID && !faceSM->DependsOn( subID )) faceID = 0; } - if ( !faceID && !isBoundary ) - continue; + // if ( !faceID && !isBoundary ) + // continue; } + if ( !faceID && !isBoundary ) + continue; } + // orient a new face according to supporting FACE orientation in shape_to_mesh - if ( !solid->IsOutsideOriented( faceID )) + if ( !isBoundary && !solid->IsOutsideOriented( faceID )) { if ( existFace ) editor.Reorient( existFace ); @@ -4865,14 +5603,15 @@ namespace } } - // split a polygon that will be used by other 3D algorithm if ( faceID && nbFaceNodes > 4 && !_grid->IsInternal( faceID ) && !_grid->IsShared( faceID ) && !_grid->IsBoundaryFace( faceID )) { - splitPolygon( newFace, vTool, iFacet, faceID, solidID, - face, editor, i+1 < bndFacets.size() ); + // split a polygon that will be used by other 3D algorithm + if ( !existFace ) + splitPolygon( newFace, vTool, iFacet, faceID, solidID, + face, editor, i+1 < bndFacets.size() ); } else { @@ -4939,6 +5678,8 @@ namespace for ( size_t i = 1; i < nodes.size(); i++ ) { + if ( mesh->FindEdge( nodes[i-1], nodes[i] )) + continue; SMDS_MeshElement* segment = mesh->AddEdge( nodes[i-1], nodes[i] ); mesh->SetMeshElementOnShape( segment, e2ff->first ); } @@ -4950,7 +5691,7 @@ namespace /*! * \brief Return created volumes and volumes that can have free facet because of * skipped small volume. Also create mesh faces on free facets - * of adjacent not-cut volumes id the result volume is too small. + * of adjacent not-cut volumes if the result volume is too small. */ void Hexahedron::getBoundaryElems( vector< const SMDS_MeshElement* > & boundaryElems ) { @@ -4997,8 +5738,9 @@ namespace if ( !faceID ) break; if ( _grid->IsInternal( faceID ) || - _grid->IsShared( faceID ) || - _grid->IsBoundaryFace( faceID )) + _grid->IsShared( faceID ) //|| + //_grid->IsBoundaryFace( faceID ) -- commented for #19887 + ) break; // create only if a new face will be used by other 3D algo } @@ -5015,7 +5757,9 @@ namespace // return created volumes for ( _volumeDef* volDef = &_volumeDefs; volDef; volDef = volDef->_next ) { - if ( volDef->_volume && !volDef->_volume->isMarked() ) + if ( volDef ->_volume && + !volDef->_volume->IsNull() && + !volDef->_volume->isMarked() ) { volDef->_volume->setIsMarked( true ); boundaryElems.push_back( volDef->_volume ); @@ -5027,6 +5771,400 @@ namespace } } + //================================================================================ + /*! + * \brief Remove edges and nodes dividing a hexa side in the case if an adjacent + * volume also sharing the dividing edge is missing due to its small side. + * Issue #19887. + */ + //================================================================================ + + void Hexahedron::removeExcessSideDivision(const vector< Hexahedron* >& allHexa) + { + if ( ! _volumeDefs.IsPolyhedron() ) + return; // not a polyhedron + + // look for a divided side adjacent to a small hexahedron + + int di[6] = { 0, 0, 0, 0,-1, 1 }; + int dj[6] = { 0, 0,-1, 1, 0, 0 }; + int dk[6] = {-1, 1, 0, 0, 0, 0 }; + + for ( int iF = 0; iF < 6; ++iF ) // loop on 6 sides of a hexahedron + { + size_t neighborIndex = _grid->CellIndex( _i + di[iF], + _j + dj[iF], + _k + dk[iF] ); + if ( neighborIndex >= allHexa.size() || + !allHexa[ neighborIndex ] || + !allHexa[ neighborIndex ]->_hasTooSmall ) + continue; + + // check if a side is divided into several polygons + for ( _volumeDef* volDef = &_volumeDefs; volDef; volDef = volDef->_next ) + { + int nbPolygons = 0, nbNodes = 0; + for ( size_t i = 0; i < volDef->_names.size(); ++i ) + if ( volDef->_names[ i ] == _hexQuads[ iF ]._name ) + { + ++nbPolygons; + nbNodes += volDef->_quantities[ i ]; + } + if ( nbPolygons < 2 ) + continue; + + // construct loops from polygons + typedef _volumeDef::_linkDef TLinkDef; + std::vector< TLinkDef* > loops; + std::vector< TLinkDef > links( nbNodes ); + for ( size_t i = 0, iN = 0, iLoop = 0; iLoop < volDef->_quantities.size(); ++iLoop ) + { + size_t nbLinks = volDef->_quantities[ iLoop ]; + if ( volDef->_names[ iLoop ] != _hexQuads[ iF ]._name ) + { + iN += nbLinks; + continue; + } + loops.push_back( & links[i] ); + for ( size_t n = 0; n < nbLinks-1; ++n, ++i, ++iN ) + { + links[i].init( volDef->_nodes[iN], volDef->_nodes[iN+1], iLoop ); + links[i].setNext( &links[i+1] ); + } + links[i].init( volDef->_nodes[iN], volDef->_nodes[iN-nbLinks+1], iLoop ); + links[i].setNext( &links[i-nbLinks+1] ); + ++i; ++iN; + } + + // look for equal links in different loops and join such loops + bool loopsJoined = false; + std::set< TLinkDef > linkSet; + for ( size_t iLoop = 0; iLoop < loops.size(); ++iLoop ) + { + TLinkDef* beg = 0; + for ( TLinkDef* l = loops[ iLoop ]; l != beg; l = l->_next ) // walk around the iLoop + { + std::pair< std::set< TLinkDef >::iterator, bool > it2new = linkSet.insert( *l ); + if ( !it2new.second ) // equal found, join loops + { + const TLinkDef* equal = &(*it2new.first); + if ( equal->_loopIndex == l->_loopIndex ) + continue; // error? + + loopsJoined = true; + + for ( size_t i = iLoop - 1; i < loops.size(); --i ) + if ( loops[ i ] && loops[ i ]->_loopIndex == equal->_loopIndex ) + loops[ i ] = 0; + + // exclude l and equal and join two loops + if ( l->_prev != equal ) + l->_prev->setNext( equal->_next ); + if ( equal->_prev != l ) + equal->_prev->setNext( l->_next ); + + if ( volDef->_quantities[ l->_loopIndex ] > 0 ) + volDef->_quantities[ l->_loopIndex ] *= -1; + if ( volDef->_quantities[ equal->_loopIndex ] > 0 ) + volDef->_quantities[ equal->_loopIndex ] *= -1; + + if ( loops[ iLoop ] == l ) + loops[ iLoop ] = l->_prev->_next; + } + beg = loops[ iLoop ]; + } + } + // update volDef + if ( loopsJoined ) + { + // set unchanged polygons + std::vector< int > newQuantities; + std::vector< _volumeDef::_nodeDef > newNodes; + vector< SMESH_Block::TShapeID > newNames; + newQuantities.reserve( volDef->_quantities.size() ); + newNodes.reserve ( volDef->_nodes.size() ); + newNames.reserve ( volDef->_names.size() ); + for ( size_t i = 0, iLoop = 0; iLoop < volDef->_quantities.size(); ++iLoop ) + { + if ( volDef->_quantities[ iLoop ] < 0 ) + { + i -= volDef->_quantities[ iLoop ]; + continue; + } + newQuantities.push_back( volDef->_quantities[ iLoop ]); + newNodes.insert( newNodes.end(), + volDef->_nodes.begin() + i, + volDef->_nodes.begin() + i + newQuantities.back() ); + newNames.push_back( volDef->_names[ iLoop ]); + i += volDef->_quantities[ iLoop ]; + } + + // set joined loops + for ( size_t iLoop = 0; iLoop < loops.size(); ++iLoop ) + { + if ( !loops[ iLoop ] ) + continue; + newQuantities.push_back( 0 ); + TLinkDef* beg = 0; + for ( TLinkDef* l = loops[ iLoop ]; l != beg; l = l->_next, ++newQuantities.back() ) + { + newNodes.push_back( l->_node1 ); + beg = loops[ iLoop ]; + } + newNames.push_back( _hexQuads[ iF ]._name ); + } + volDef->_quantities.swap( newQuantities ); + volDef->_nodes.swap( newNodes ); + volDef->_names.swap( newNames ); + } + } // loop on volDef's + } // loop on hex sides + + return; + } // removeExcessSideDivision() + + + //================================================================================ + /*! + * \brief Remove nodes splitting Cartesian cell edges in the case if a node + * is used in every cells only by two polygons sharing the edge + * Issue #19887. + */ + //================================================================================ + + void Hexahedron::removeExcessNodes(vector< Hexahedron* >& allHexa) + { + if ( ! _volumeDefs.IsPolyhedron() ) + return; // not a polyhedron + + typedef vector< _volumeDef::_nodeDef >::iterator TNodeIt; + vector< int > nodesInPoly[ 4 ]; // node index in _volumeDefs._nodes + vector< int > volDefInd [ 4 ]; // index of a _volumeDefs + Hexahedron* hexa [ 4 ]; + int i,j,k, cellIndex, iLink = 0, iCellLink; + for ( int iDir = 0; iDir < 3; ++iDir ) + { + CellsAroundLink fourCells( _grid, iDir ); + for ( int iL = 0; iL < 4; ++iL, ++iLink ) // 4 links in a direction + { + _Link& link = _hexLinks[ iLink ]; + fourCells.Init( _i, _j, _k, iLink ); + + for ( size_t iP = 0; iP < link._fIntPoints.size(); ++iP ) // loop on nodes on the link + { + bool nodeRemoved = true; + _volumeDef::_nodeDef node; node._intPoint = link._fIntPoints[iP]; + + for ( size_t i = 0, nb = _volumeDefs.size(); i < nb && nodeRemoved; ++i ) + if ( _volumeDef* vol = _volumeDefs.at( i )) + nodeRemoved = + ( std::find( vol->_nodes.begin(), vol->_nodes.end(), node ) == vol->_nodes.end() ); + if ( nodeRemoved ) + continue; // node already removed + + // check if a node encounters zero or two times in 4 cells sharing iLink + // if so, the node can be removed from the cells + bool nodeIsOnEdge = true; + int nbPolyhedraWithNode = 0; + for ( int iC = 0; iC < 4; ++iC ) // loop on 4 cells sharing a link + { + nodesInPoly[ iC ].clear(); + volDefInd [ iC ].clear(); + hexa [ iC ] = 0; + if ( !fourCells.GetCell( iC, i,j,k, cellIndex, iCellLink )) + continue; + hexa[ iC ] = allHexa[ cellIndex ]; + if ( !hexa[ iC ]) + continue; + for ( size_t i = 0, nb = hexa[ iC ]->_volumeDefs.size(); i < nb; ++i ) + if ( _volumeDef* vol = hexa[ iC ]->_volumeDefs.at( i )) + { + for ( TNodeIt nIt = vol->_nodes.begin(); nIt != vol->_nodes.end(); ++nIt ) + { + nIt = std::find( nIt, vol->_nodes.end(), node ); + if ( nIt != vol->_nodes.end() ) + { + nodesInPoly[ iC ].push_back( std::distance( vol->_nodes.begin(), nIt )); + volDefInd [ iC ].push_back( i ); + } + else + break; + } + nbPolyhedraWithNode += ( !nodesInPoly[ iC ].empty() ); + } + if ( nodesInPoly[ iC ].size() != 0 && + nodesInPoly[ iC ].size() != 2 ) + { + nodeIsOnEdge = false; + break; + } + } // loop on 4 cells + + // remove nodes from polyhedra + if ( nbPolyhedraWithNode > 0 && nodeIsOnEdge ) + { + for ( int iC = 0; iC < 4; ++iC ) // loop on 4 cells sharing the link + { + if ( nodesInPoly[ iC ].empty() ) + continue; + for ( int i = volDefInd[ iC ].size() - 1; i >= 0; --i ) + { + _volumeDef* vol = hexa[ iC ]->_volumeDefs.at( volDefInd[ iC ][ i ]); + int nIndex = nodesInPoly[ iC ][ i ]; + // decrement _quantities + for ( size_t iQ = 0; iQ < vol->_quantities.size(); ++iQ ) + if ( nIndex < vol->_quantities[ iQ ]) + { + vol->_quantities[ iQ ]--; + break; + } + else + { + nIndex -= vol->_quantities[ iQ ]; + } + vol->_nodes.erase( vol->_nodes.begin() + nodesInPoly[ iC ][ i ]); + + if ( i == 0 && + vol->_nodes.size() == 6 * 4 && + vol->_quantities.size() == 6 ) // polyhedron becomes hexahedron? + { + bool allQuads = true; + for ( size_t iQ = 0; iQ < vol->_quantities.size() && allQuads; ++iQ ) + allQuads = ( vol->_quantities[ iQ ] == 4 ); + if ( allQuads ) + { + // set side nodes as this: bottom, top, top, ... + int iTop = 0, iBot = 0; // side indices + for ( int iS = 0; iS < 6; ++iS ) + { + if ( vol->_names[ iS ] == SMESH_Block::ID_Fxy0 ) + iBot = iS; + if ( vol->_names[ iS ] == SMESH_Block::ID_Fxy1 ) + iTop = iS; + } + if ( iBot != 0 ) + { + if ( iTop == 0 ) + { + std::copy( vol->_nodes.begin(), + vol->_nodes.begin() + 4, + vol->_nodes.begin() + 4 ); + iTop = 1; + } + std::copy( vol->_nodes.begin() + 4 * iBot, + vol->_nodes.begin() + 4 * ( iBot + 1), + vol->_nodes.begin() ); + } + if ( iTop != 1 ) + std::copy( vol->_nodes.begin() + 4 * iTop, + vol->_nodes.begin() + 4 * ( iTop + 1), + vol->_nodes.begin() + 4 ); + + std::copy( vol->_nodes.begin() + 4, + vol->_nodes.begin() + 8, + vol->_nodes.begin() + 8 ); + // set up top facet nodes by comparing their uvw with bottom nodes + E_IntersectPoint ip[8]; + for ( int iN = 0; iN < 8; ++iN ) + { + SMESH_NodeXYZ p = vol->_nodes[ iN ].Node(); + _grid->ComputeUVW( p, ip[ iN ]._uvw ); + } + const double tol2 = _grid->_tol * _grid->_tol; + for ( int iN = 0; iN < 4; ++iN ) + { + gp_Pnt2d pBot( ip[ iN ]._uvw[0], ip[ iN ]._uvw[1] ); + for ( int iT = 4; iT < 8; ++iT ) + { + gp_Pnt2d pTop( ip[ iT ]._uvw[0], ip[ iT ]._uvw[1] ); + if ( pBot.SquareDistance( pTop ) < tol2 ) + { + // vol->_nodes[ iN + 4 ]._node = ip[ iT ]._node; + // vol->_nodes[ iN + 4 ]._intPoint = 0; + vol->_nodes[ iN + 4 ] = vol->_nodes[ iT + 4 ]; + break; + } + } + } + vol->_nodes.resize( 8 ); + vol->_quantities.clear(); + //vol->_names.clear(); + } + } + } // loop on _volumeDefs + } // loop on 4 cell abound a link + } // if ( nodeIsOnEdge ) + } // loop on intersection points of a link + } // loop on 4 links of a direction + } // loop on 3 directions + + return; + + } // removeExcessNodes() + + //================================================================================ + /*! + * \brief [Issue #19913] Modify _hexLinks._splits to prevent creating overlapping volumes + */ + //================================================================================ + + void Hexahedron::preventVolumesOverlapping() + { + // Cut off a quadrangle corner if two links sharing the corner + // are shared by same two solids, in this case each of solids gets + // a triangle for it-self. + std::vector< TGeomID > soIDs[4]; + for ( int iF = 0; iF < 6; ++iF ) // loop on 6 sides of a hexahedron + { + _Face& quad = _hexQuads[ iF ] ; + + int iFOpposite = iF + ( iF % 2 ? -1 : 1 ); + _Face& quadOpp = _hexQuads[ iFOpposite ] ; + + int nbSides = 0, nbSidesOpp = 0; + for ( int iE = 0; iE < 4; ++iE ) // loop on 4 sides of a quadrangle + { + nbSides += ( quad._links [ iE ].NbResultLinks() > 0 ); + nbSidesOpp += ( quadOpp._links[ iE ].NbResultLinks() > 0 ); + } + if ( nbSides < 4 || nbSidesOpp != 2 ) + continue; + + for ( int iE = 0; iE < 4; ++iE ) + { + soIDs[ iE ].clear(); + _Node* n = quad._links[ iE ].FirstNode(); + if ( n->_intPoint && n->_intPoint->_faceIDs.size() ) + soIDs[ iE ] = _grid->GetSolidIDs( n->_intPoint->_faceIDs[0] ); + } + if ((( soIDs[0].size() >= 2 ) + + ( soIDs[1].size() >= 2 ) + + ( soIDs[2].size() >= 2 ) + + ( soIDs[3].size() >= 2 ) ) < 3 ) + continue; + + bool done = false; + for ( int i = 0; i < 4; ++i ) + { + int i1 = _grid->_helper->WrapIndex( i + 1, 4 ); + int i2 = _grid->_helper->WrapIndex( i + 2, 4 ); + int i3 = _grid->_helper->WrapIndex( i + 3, 4 ); + if ( soIDs[i1].size() == 2 && soIDs[i ] != soIDs[i1] && + soIDs[i2].size() == 2 && soIDs[i1] == soIDs[i2] && + soIDs[i3].size() == 2 && soIDs[i2] == soIDs[i3] ) + { + quad._links[ i1 ]._link->_splits.clear(); + quad._links[ i2 ]._link->_splits.clear(); + done = true; + break; + } + } + if ( done ) + break; + } + return; + } // preventVolumesOverlapping() + //================================================================================ /*! * \brief Set to _hexLinks a next portion of splits located on one side of INTERNAL FACEs @@ -5357,7 +6495,10 @@ bool StdMeshers_Cartesian_3D::Compute(SMESH_Mesh & theMesh, { multiset< F_IntersectPoint >::iterator ip = lines[i]._intPoints.begin(); for ( ; ip != lines[i]._intPoints.end(); ++ip ) - if ( ip->_node && ip->_node->NbInverseElements() == 0 && !ip->_node->isMarked() ) + if ( ip->_node && + !ip->_node->IsNull() && + ip->_node->NbInverseElements() == 0 && + !ip->_node->isMarked() ) { nodesToRemove.push_back( ip->_node ); ip->_node->setIsMarked( true ); @@ -5366,7 +6507,9 @@ bool StdMeshers_Cartesian_3D::Compute(SMESH_Mesh & theMesh, } // get grid nodes for ( size_t i = 0; i < grid._nodes.size(); ++i ) - if ( grid._nodes[i] && grid._nodes[i]->NbInverseElements() == 0 && + if ( grid._nodes[i] && + !grid._nodes[i]->IsNull() && + grid._nodes[i]->NbInverseElements() == 0 && !grid._nodes[i]->isMarked() ) { nodesToRemove.push_back( grid._nodes[i] ); @@ -5395,9 +6538,9 @@ bool StdMeshers_Cartesian_3D::Compute(SMESH_Mesh & theMesh, */ //============================================================================= -bool StdMeshers_Cartesian_3D::Evaluate(SMESH_Mesh & theMesh, - const TopoDS_Shape & theShape, - MapShapeNbElems& theResMap) +bool StdMeshers_Cartesian_3D::Evaluate(SMESH_Mesh & /*theMesh*/, + const TopoDS_Shape & /*theShape*/, + MapShapeNbElems& /*theResMap*/) { // TODO // std::vector aResVec(SMDSEntity_Last); @@ -5451,11 +6594,11 @@ namespace // -------------------------------------------------------------------------------- // unsetting _alwaysComputed flag if "Cartesian_3D" was removed // - virtual void ProcessEvent(const int event, + virtual void ProcessEvent(const int /*event*/, const int eventType, SMESH_subMesh* subMeshOfSolid, - SMESH_subMeshEventListenerData* data, - const SMESH_Hypothesis* hyp = 0) + SMESH_subMeshEventListenerData* /*data*/, + const SMESH_Hypothesis* /*hyp*/ = 0) { if ( eventType == SMESH_subMesh::COMPUTE_EVENT ) {