From 36b825e9340bdf7c581df556948ff1903d1eac5c Mon Sep 17 00:00:00 2001 From: Konstantin Leontev Date: Tue, 30 Jul 2024 17:36:57 +0100 Subject: [PATCH] [bos #42217][EDF 28921] Horseshoe with bodyfitting. Nodes linked to a null nodes without any possible splits now being cleared befor splits initialization. --- src/StdMeshers/StdMeshers_Cartesian_3D.cxx | 232 +++++++++++++++++++-- 1 file changed, 211 insertions(+), 21 deletions(-) diff --git a/src/StdMeshers/StdMeshers_Cartesian_3D.cxx b/src/StdMeshers/StdMeshers_Cartesian_3D.cxx index 4a0f7070e..34bf2fb61 100644 --- a/src/StdMeshers/StdMeshers_Cartesian_3D.cxx +++ b/src/StdMeshers/StdMeshers_Cartesian_3D.cxx @@ -203,6 +203,8 @@ namespace const TGeomID theUndefID = 1e+9; + const std::string debugSepLine = "\n===========================================\n"; + //============================================================================= // Definitions of internal utils // -------------------------------------------------------------------------- @@ -774,12 +776,23 @@ namespace _node = _intPoint->_node = node; } + void clear() + { + _node = nullptr; + _boundaryCornerNode = nullptr; + _intPoint = nullptr; + _usedInFace = nullptr; + _isInternalFlags = IS_NOT_INTERNAL; + } + friend std::ostream& operator<<(std::ostream& os, const _Node& node) { if (node._node) - node._node->Print(os); + node._node->Print(os); // mesh node at hexahedron corner + else if (node._intPoint && node._intPoint->_node) + node._intPoint->_node->Print(os); // intersection point else - os << "mesh node at hexahedron corner is null" << '\n'; + os << "mesh node is null" << '\n'; return os; } @@ -787,18 +800,33 @@ namespace // -------------------------------------------------------------------------------- struct _Link // link connecting two _Node's { - _Node* _nodes[2]; - _Face* _faces[2]; // polygons sharing a link + static const std::size_t arrSize = 2; + + _Node* _nodes[arrSize]; + _Face* _faces[arrSize]; // polygons sharing a link vector< const F_IntersectPoint* > _fIntPoints; // GridLine intersections with FACEs vector< _Node* > _fIntNodes; // _Node's at _fIntPoints vector< _Link > _splits; _Link(): _faces{ 0, 0 } {} + void clear() + { + for (std::size_t i = 0; i < arrSize; ++i) + { + if (_nodes[i]) + _nodes[i]->clear(); + } + + _fIntPoints.clear(); + _fIntNodes.clear(); + _splits.clear(); + } + friend std::ostream& operator<<(std::ostream& os, const _Link& link) { os << "Link:\n"; - for (int i = 0; i < 2; ++i) + for (std::size_t i = 0; i < arrSize; ++i) { if (link._nodes[i]) os << *link._nodes[i]; @@ -806,6 +834,10 @@ namespace os << "link node with index " << i << " is null" << '\n'; } + os << "_fIntPoints: " << link._fIntPoints.size() << '\n'; + os << "_fIntNodes: " << link._fIntNodes.size() << '\n'; + os << "_splits: " << link._splits.size() << '\n'; + return os; } }; @@ -958,6 +990,34 @@ namespace _polyLinks.push_back( l ); _links.push_back( _OrientedLink( &_polyLinks.back() )); } + + friend std::ostream& operator<<(std::ostream& os, const _Face& face) + { + os << "Face " << face._name << '\n'; + + os << "Links on GridLines: \n"; + for (const auto& link : face._links) + { + os << link; + } + + os << "Links added to close a polygonal face: \n"; + for (const auto& link : face._polyLinks) + { + os << link; + } + + os << "Nodes at intersection with EDGEs: \n"; + for (const auto node : face._eIntNodes) + { + if (node) + { + os << *node; + } + } + + return os; + } }; // -------------------------------------------------------------------------------- struct _volumeDef // holder of nodes of a volume mesh element @@ -1037,9 +1097,12 @@ namespace }; // topology of a hexahedron - _Node _hexNodes [8]; - _Link _hexLinks [12]; - _Face _hexQuads [6]; + static const std::size_t HEX_NODES_NUM = 8; + static const std::size_t HEX_LINKS_NUM = 12; + static const std::size_t HEX_QUADS_NUM = 6; + _Node _hexNodes [HEX_NODES_NUM]; + _Link _hexLinks [HEX_LINKS_NUM]; + _Face _hexQuads [HEX_QUADS_NUM]; // faces resulted from hexahedron intersection vector< _Face > _polygons; @@ -1074,6 +1137,8 @@ namespace Hexahedron(const Hexahedron& other, size_t i, size_t j, size_t k, int cellID ); void init( size_t i, size_t j, size_t k, const Solid* solid=0 ); void init( size_t i ); + void clearNodesLinkedToNull(const Solid* solid, SMESH_MesherHelper& helper); + bool isSplittedLink(const Solid* solid, SMESH_MesherHelper& helper, const Hexahedron::_Link& linkIn) const; void setIJK( size_t i ); bool compute( const Solid* solid, const IsInternalFlag intFlag ); size_t getSolids( TGeomID ids[] ); @@ -2726,12 +2791,108 @@ namespace init( _i, _j, _k ); } + //================================================================================ + /*! + * \brief Clear nodes linked to null by not splitted links + */ + void Hexahedron::clearNodesLinkedToNull(const Solid* solid, SMESH_MesherHelper& helper) + { + for (std::size_t i = 0; i < HEX_LINKS_NUM; ++i) + { + _Link& link = _hexLinks[i]; + + if (isSplittedLink(solid, helper, link)) + continue; + + // Links where both nodes are valid should have itself as a split. + // Check if we have at least one null node here for a chance if we have + // some unexpected edge case. + _Node* n1 = link._nodes[0]; + _Node* n2 = link._nodes[1]; + + assert(!n1->_node || !n2->_node); + if (!n1->_node || !n2->_node) + { + link.clear(); + } + } + } + + //================================================================================ + /*! + * \brief Returns true if the link will be splitted on the Hexaedron initialization + */ + bool Hexahedron::isSplittedLink(const Solid* solid, SMESH_MesherHelper& helper, const Hexahedron::_Link& linkIn) const + { + _Link split; + std::vector<_Node> intNodes; + + _Link link = linkIn; + link._fIntNodes.clear(); + link._fIntNodes.reserve( link._fIntPoints.size() ); + + for ( size_t i = 0; i < link._fIntPoints.size(); ++i ) + if ( solid->ContainsAny( link._fIntPoints[i]->_faceIDs )) + { + intNodes.push_back( _Node( 0, link._fIntPoints[i] )); + link._fIntNodes.push_back( & intNodes.back() ); + } + + link._splits.clear(); + split._nodes[ 0 ] = link._nodes[0]; + bool isOut = ( ! link._nodes[0]->Node() ); + bool checkTransition; + for ( size_t i = 0; i < link._fIntNodes.size(); ++i ) + { + const bool isGridNode = ( ! link._fIntNodes[i]->Node() ); + if ( !isGridNode ) // intersection non-coincident with a grid node + { + if ( split._nodes[ 0 ]->Node() && !isOut ) + { + return true; + } + split._nodes[ 0 ] = link._fIntNodes[i]; + checkTransition = true; + } + else // FACE intersection coincident with a grid node (at link ends) + { + checkTransition = ( i == 0 && link._nodes[0]->Node() ); + } + if ( checkTransition ) + { + const vector< TGeomID >& faceIDs = link._fIntNodes[i]->_intPoint->_faceIDs; + if ( _grid->IsInternal( faceIDs.back() )) + isOut = false; + else if ( faceIDs.size() > 1 || _eIntPoints.size() > 0 ) + isOut = isOutPoint( link, i, helper, solid ); + else + { + bool okTransi = _grid->IsCorrectTransition( faceIDs[0], solid ); + switch ( link._fIntNodes[i]->FaceIntPnt()->_transition ) { + case Trans_OUT: isOut = okTransi; break; + case Trans_IN : isOut = !okTransi; break; + default: + isOut = isOutPoint( link, i, helper, solid ); + } + } + } + } + if ( link._nodes[ 1 ]->Node() && split._nodes[ 0 ]->Node() && !isOut ) + { + return true; + } + + return false; + } + //================================================================================ /*! * \brief Initializes its data by given grid cell nodes and intersections */ void Hexahedron::init( size_t i, size_t j, size_t k, const Solid* solid ) { + MESSAGE(debugSepLine << "START Hexahedron::init()" << debugSepLine); + _i = i; _j = j; _k = k; bool isCompute = solid; @@ -2780,6 +2941,8 @@ namespace // this method can be called in parallel, so use own helper SMESH_MesherHelper helper( *_grid->_helper->GetMesh() ); + clearNodesLinkedToNull(solid, helper); + // Create sub-links (_Link::_splits) by splitting links with _Link::_fIntPoints // --------------------------------------------------------------- _Link split; @@ -2968,8 +3131,8 @@ namespace } } } - return; + MESSAGE(debugSepLine << "END Hexahedron::init()" << debugSepLine); } // init( _i, _j, _k ) //================================================================================ @@ -3084,23 +3247,26 @@ namespace */ bool Hexahedron::collectSplits(std::vector<_OrientedLink>& splits, const _Face& quad, _Face* polygon, int quadIndex) { - MESSAGE("Collect splits..."); + MESSAGE(debugSepLine << "START Collect splits..." << debugSepLine); splits.clear(); for ( int iE = 0; iE < 4; ++iE ) // loop on 4 sides of a quadrangle for ( size_t iS = 0; iS < quad._links[ iE ].NbResultLinks(); ++iS ) + { splits.push_back( quad._links[ iE ].ResultLink( iS )); + // SCRUTE(splits.back()); + } if ( splits.size() == 4 && isQuadOnFace( quadIndex )) // check if a quad on FACE is not split { - MESSAGE("The quad on FACE is not split. Swap splits with polygon links."); + MESSAGE(debugSepLine << "END. The quad on FACE is not split. Swap splits with polygon links." << debugSepLine); polygon->_links.swap( splits ); return false; } - MESSAGE("Collected splits: " << splits.size()); + MESSAGE(debugSepLine << "Collected splits: " << splits.size() << debugSepLine); return true; } @@ -3132,7 +3298,7 @@ namespace void Hexahedron::connectPolygonLinks( const Solid* solid, _Face* polygon, _Face& quad, vector<_Node*>& chainNodes, std::vector<_OrientedLink>& splits, const bool toCheckSideDivision) { - MESSAGE("Connect polygon links..."); + MESSAGE(debugSepLine << "START connectPolygonLinks()..." << debugSepLine); const std::set concaveFaces = getConcaveFaces(solid); // to avoid connecting nodes laying on them @@ -3164,6 +3330,7 @@ void Hexahedron::connectPolygonLinks( polygon = createPolygon(quad._name); } polygon->_links.push_back( splits[ iS ] ); + MESSAGE("Split link was added to a polygon: " << splits[iS]); splits[ iS++ ]._link = 0; --nbSplits; @@ -3200,10 +3367,11 @@ void Hexahedron::connectPolygonLinks( else if ( n1 != n2 ) { // try to connect to intersections with EDGEs - MESSAGE("n1 != n2 -> try to connect to intersections with EDGEs"); + MESSAGE("n1 != n2..."); if ( quad._eIntNodes.size() > nbUsedEdgeNodes && findChain( n2, n1, quad, chainNodes )) { + MESSAGE("try to connect to intersections with EDGEs"); for ( size_t i = 1; i < chainNodes.size(); ++i ) { polygon->AddPolyLink( chainNodes[i-1], chainNodes[i] ); @@ -3237,6 +3405,7 @@ void Hexahedron::connectPolygonLinks( { if ( n2 != foundSplit.FirstNode() ) { + MESSAGE("Add link n2 -> foundSplit.FirstNode()"); polygon->AddPolyLink( n2, foundSplit.FirstNode() ); n2 = foundSplit.FirstNode(); } @@ -3271,7 +3440,7 @@ void Hexahedron::connectPolygonLinks( } } // while ( nbSplits > 0 ) - MESSAGE("Polygon's links were connected"); + MESSAGE(debugSepLine << "END connectPolygonLinks(): " << *polygon << debugSepLine); } //================================================================================ @@ -3280,6 +3449,8 @@ void Hexahedron::connectPolygonLinks( */ std::vector Hexahedron::getChainNodes(const Solid* solid, const IsInternalFlag intFlag) { + MESSAGE(debugSepLine << "START Get chain nodes..." << debugSepLine); + std::vector< _OrientedLink > splits; std::vector<_Node*> chainNodes; @@ -3289,6 +3460,7 @@ void Hexahedron::connectPolygonLinks( { _Face& quad = _hexQuads[ iF ] ; _Face* polygon = createPolygon(quad._name); + SCRUTE(quad); if (!collectSplits(splits, quad, polygon, iF)) continue; // goto the next quad @@ -3300,6 +3472,7 @@ void Hexahedron::connectPolygonLinks( } } // loop on 6 hexahedron sides + MESSAGE(debugSepLine << "END Collected " << chainNodes.size() << " chain nodes" << debugSepLine); return chainNodes; } @@ -3329,14 +3502,17 @@ void Hexahedron::connectPolygonLinks( */ void Hexahedron::addPolygonsToLinks() { + MESSAGE(debugSepLine << "START addPolygonsToLinks()..." << debugSepLine); for (_Face& polygon : _polygons) { + // MESSAGE("START add polygon to its links : " << polygon << debugSepLine); for (_OrientedLink& link : polygon._links) { link.AddFace( &polygon ); link.FirstNode()->_usedInFace = &polygon; } } + MESSAGE(debugSepLine << "END addPolygonsToLinks()" << debugSepLine); } //================================================================================ @@ -3393,7 +3569,7 @@ void Hexahedron::connectPolygonLinks( */ bool Hexahedron::createPolygons(const bool hasEdgeIntersections, const IsInternalFlag intFlag) { - MESSAGE("Create polygons from free links..."); + MESSAGE(debugSepLine << "START createPolygons()..." << debugSepLine); // find free links vector< _OrientedLink* > freeLinks = getFreeLinks(); @@ -3408,6 +3584,7 @@ void Hexahedron::connectPolygonLinks( std::vector< TGeomID > faces; TGeomID curFace = 0; const size_t nbQuadPolygons = _polygons.size(); + SCRUTE(nbQuadPolygons); E_IntersectPoint ipTmp; std::map< TGeomID, std::vector< const B_IntersectPoint* > > tmpAddedFace; // face added to _intPoint @@ -3423,6 +3600,7 @@ void Hexahedron::connectPolygonLinks( _polygons[ iPolygon ]._links.reserve( 20 ); } _Face& polygon = _polygons[ iPolygon ]; + SCRUTE(polygon); _OrientedLink* curLink = 0; _Node* curNode; @@ -3450,6 +3628,7 @@ void Hexahedron::connectPolygonLinks( polygon._links.push_back( *curLink ); } } while ( curLink ); + SCRUTE(polygon); } else // there are intersections with EDGEs { @@ -3605,6 +3784,7 @@ void Hexahedron::connectPolygonLinks( { polygon._links[ iL ].AddFace( &polygon ); polygon._links[ iL ].Reverse(); + SCRUTE(polygon); } if ( /*hasEdgeIntersections &&*/ iPolygon == _polygons.size() - 1 ) { @@ -3729,6 +3909,7 @@ void Hexahedron::connectPolygonLinks( if ( _hexNodes[ i ]._intPoint == &ipTmp ) _hexNodes[ i ]._intPoint = 0; + MESSAGE(debugSepLine << "END createPolygons()" << debugSepLine); return !_hasTooSmall; // too small volume } @@ -3777,7 +3958,7 @@ void Hexahedron::connectPolygonLinks( */ bool Hexahedron::createVolume(const Solid* solid) { - MESSAGE("Create a volume..."); + MESSAGE(debugSepLine << "START createVolume()" << debugSepLine); _volumeDefs._nodes.clear(); _volumeDefs._quantities.clear(); @@ -3785,9 +3966,12 @@ bool Hexahedron::createVolume(const Solid* solid) // create a classic cell if possible int nbPolygons = 0; + SCRUTE(_polygons.size()); for ( size_t iF = 0; iF < _polygons.size(); ++iF ) nbPolygons += (_polygons[ iF ]._links.size() > 2 ); + SCRUTE(nbPolygons); + //const int nbNodes = _nbCornerNodes + nbIntersections; int nbNodes = 0; for ( size_t i = 0; i < 8; ++i ) @@ -3795,6 +3979,7 @@ bool Hexahedron::createVolume(const Solid* solid) for ( size_t i = 0; i < _intNodes.size(); ++i ) nbNodes += _intNodes[ i ].IsUsedInFace(); + SCRUTE(nbNodes); bool isClassicElem = false; if ( nbNodes == 8 && nbPolygons == 6 ) isClassicElem = addHexa(); else if ( nbNodes == 4 && nbPolygons == 4 ) isClassicElem = addTetra(); @@ -3807,17 +3992,18 @@ bool Hexahedron::createVolume(const Solid* solid) for ( size_t iF = 0; iF < _polygons.size(); ++iF ) { const size_t nbLinks = _polygons[ iF ]._links.size(); - SCRUTE(nbLinks); if ( nbLinks < 3 ) continue; _volumeDefs._quantities.push_back( nbLinks ); _volumeDefs._names.push_back( _polygons[ iF ]._name ); - SCRUTE(_polygons[ iF ]._name); for ( size_t iL = 0; iL < nbLinks; ++iL ) _volumeDefs._nodes.push_back( _polygons[ iF ]._links[ iL ].FirstNode() ); + + // SCRUTE(_polygons[ iF ]); } } _volumeDefs._solidID = solid->ID(); + MESSAGE(debugSepLine << "END createVolume()" << debugSepLine); return !_volumeDefs._nodes.empty(); } @@ -3826,7 +4012,7 @@ bool Hexahedron::createVolume(const Solid* solid) */ bool Hexahedron::compute( const Solid* solid, const IsInternalFlag intFlag ) { - MESSAGE("Compute solid " << solid->ID() << " with internal flag " << intFlag << " ..."); + MESSAGE(debugSepLine << "START Compute solid " << solid->ID() << " ..." << debugSepLine); // Reset polygons _polygons.clear(); @@ -3858,7 +4044,11 @@ bool Hexahedron::createVolume(const Solid* solid) // Fix polygons' names setNamesForNoNamePolygons(); - return createVolume(solid); + const bool isVolumeCreated = createVolume(solid); + SCRUTE(isVolumeCreated); + MESSAGE(debugSepLine << "END Compute solid " << solid->ID() << debugSepLine); + + return isVolumeCreated; } template -- 2.39.2