X-Git-Url: http://git.salome-platform.org/gitweb/?p=modules%2Fsmesh.git;a=blobdiff_plain;f=src%2FStdMeshers%2FStdMeshers_RadialQuadrangle_1D2D.cxx;h=ff5906d53344e51b150f7562ce4135e70ffce33b;hp=8d165c8e5f836ed11efaeb23afc630987a1a860e;hb=499f29d24922cec66e41b41a0039a954993bc6df;hpb=776a610c849651350f344157ec7fbb170ca3e8cc diff --git a/src/StdMeshers/StdMeshers_RadialQuadrangle_1D2D.cxx b/src/StdMeshers/StdMeshers_RadialQuadrangle_1D2D.cxx index 8d165c8e5..ff5906d53 100644 --- a/src/StdMeshers/StdMeshers_RadialQuadrangle_1D2D.cxx +++ b/src/StdMeshers/StdMeshers_RadialQuadrangle_1D2D.cxx @@ -1,4 +1,4 @@ -// Copyright (C) 2007-2016 CEA/DEN, EDF R&D, OPEN CASCADE +// Copyright (C) 2007-2022 CEA/DEN, EDF R&D, OPEN CASCADE // // This library is free software; you can redistribute it and/or // modify it under the terms of the GNU Lesser General Public @@ -74,9 +74,8 @@ using namespace std; //======================================================================= StdMeshers_RadialQuadrangle_1D2D::StdMeshers_RadialQuadrangle_1D2D(int hypId, - int studyId, SMESH_Gen* gen) - :StdMeshers_Quadrangle_2D( hypId, studyId, gen ) + :StdMeshers_Quadrangle_2D( hypId, gen ) { _name = "RadialQuadrangle_1D2D"; _shapeType = (1 << TopAbs_FACE); // 1 bit per shape type @@ -164,38 +163,45 @@ namespace static TEdgeMarker theEdgeMarker; return &theEdgeMarker; } - //! Clear face sumbesh if something happens on edges - void ProcessEvent(const int event, + //! Clear edge sumbesh if something happens on face + void ProcessEvent(const int /*event*/, const int eventType, - SMESH_subMesh* edgeSubMesh, - EventListenerData* data, + SMESH_subMesh* /*faceSubMesh*/, + EventListenerData* edgesHolder, const SMESH_Hypothesis* /*hyp*/) { - if ( data && !data->mySubMeshes.empty() && eventType == SMESH_subMesh::ALGO_EVENT) + if ( edgesHolder && eventType == SMESH_subMesh::ALGO_EVENT) { - ASSERT( data->mySubMeshes.front() != edgeSubMesh ); - SMESH_subMesh* faceSubMesh = data->mySubMeshes.front(); - faceSubMesh->ComputeStateEngine( SMESH_subMesh::CLEAN ); + std::list::iterator smIt = edgesHolder->mySubMeshes.begin(); + for ( ; smIt != edgesHolder->mySubMeshes.end(); ++smIt ) + { + SMESH_subMesh* edgeSM = *smIt; + edgeSM->ComputeStateEngine( SMESH_subMesh::CLEAN ); + } } } - }; - - //================================================================================ - /*! - * \brief Mark an edge as computed by StdMeshers_RadialQuadrangle_1D2D - */ - //================================================================================ - - void markEdgeAsComputedByMe(const TopoDS_Edge& edge, SMESH_subMesh* faceSubMesh) - { - if ( SMESH_subMesh* edgeSM = faceSubMesh->GetFather()->GetSubMeshContaining( edge )) + //! Store edge SMESH_subMesh'es computed by the algo + static void markEdge( const TopoDS_Edge& edge, SMESH_subMesh* faceSM ) { - if ( !edgeSM->GetEventListenerData( TEdgeMarker::getListener() )) - faceSubMesh->SetEventListener( TEdgeMarker::getListener(), - SMESH_subMeshEventListenerData::MakeData(faceSubMesh), - edgeSM); + if ( SMESH_subMesh* edgeSM = faceSM->GetFather()->GetSubMeshContaining( edge )) + { + EventListenerData* edgesHolder = faceSM->GetEventListenerData( getListener() ); + if ( edgesHolder ) + { + std::list::iterator smIt = std::find( edgesHolder->mySubMeshes.begin(), + edgesHolder->mySubMeshes.end(), + edgeSM ); + if ( smIt == edgesHolder->mySubMeshes.end() ) + edgesHolder->mySubMeshes.push_back( edgeSM ); + } + else + { + edgesHolder = SMESH_subMeshEventListenerData::MakeData( edgeSM ); + faceSM->SetEventListener( TEdgeMarker::getListener(), edgesHolder, faceSM ); + } + } } - } + }; //================================================================================ /*! @@ -220,12 +226,19 @@ namespace if ( nbWire > 2 || nbEdgesInWire.front() < 1 ) return 0; // remove degenerated EDGEs + TopTools_MapOfShape degenVV; list::iterator edge = edges.begin(); while ( edge != edges.end() ) if ( SMESH_Algo::isDegenerated( *edge )) + { + degenVV.Add( SMESH_MesherHelper::IthVertex( 0, *edge )); + degenVV.Add( SMESH_MesherHelper::IthVertex( 1, *edge )); edge = edges.erase( edge ); + } else + { ++edge; + } int nbEdges = edges.size(); // find VERTEXes between continues EDGEs @@ -321,9 +334,25 @@ namespace double len1 = (++l2i)->first; double len2 = (++l2i)->first; if ( len1 - len0 > len2 - len1 ) - deviation2sideInd.insert( make_pair( 0., len2sideInd.begin()->second )); + deviation2sideInd.insert( std::make_pair( 0., len2sideInd.begin()->second )); else - deviation2sideInd.insert( make_pair( 0., len2sideInd.rbegin()->second )); + deviation2sideInd.insert( std::make_pair( 0., len2sideInd.rbegin()->second )); + } + + double minDevi = deviation2sideInd.begin()->first; + int iMinCurv = deviation2sideInd.begin()->second; + if ( sides.size() == 3 && degenVV.Size() == 1 && + minDevi / sides[ iMinCurv ]->Length() > 1e-3 ) + { + // a triangle with curved sides and a degenerated EDGE (IPAL54585); + // use a side opposite to the degenerated EDGE as an elliptic one + for ( size_t iS = 0; iS < sides.size(); ++iS ) + if ( degenVV.Contains( sides[ iS ]->FirstVertex() )) + { + deviation2sideInd.clear(); + deviation2sideInd.insert( std::make_pair( 0.,( iS + 1 ) % sides.size() )); + break; + } } int iCirc = deviation2sideInd.rbegin()->second; @@ -373,9 +402,9 @@ namespace */ //================================================================================ - bool isCornerInsideCircle(const StdMeshers_FaceSidePtr& CircSide, - const StdMeshers_FaceSidePtr& LinSide1, - const StdMeshers_FaceSidePtr& LinSide2) + bool isCornerInsideCircle(const StdMeshers_FaceSidePtr& /*CircSide*/, + const StdMeshers_FaceSidePtr& /*LinSide1*/, + const StdMeshers_FaceSidePtr& /*LinSide2*/) { // if ( CircSide && LinSide1 && LinSide2 ) // { @@ -410,7 +439,7 @@ namespace { // find the center and a point most distant from it - double maxDist = 0, normPar; + double maxDist = 0, normPar = 0; gp_XY uv1, uv2; for ( int i = 0; i < 32; ++i ) { @@ -631,7 +660,7 @@ public: const int myID = -1001; TNodeDistributor* myHyp = dynamic_cast( aMesh.GetHypothesis( myID )); if ( !myHyp ) - myHyp = new TNodeDistributor( myID, 0, aMesh.GetGen() ); + myHyp = new TNodeDistributor( myID, aMesh.GetGen() ); return myHyp; } // ----------------------------------------------------------------------------- @@ -726,8 +755,8 @@ public: } protected: // ----------------------------------------------------------------------------- - TNodeDistributor( int hypId, int studyId, SMESH_Gen* gen) - : StdMeshers_Regular_1D( hypId, studyId, gen) + TNodeDistributor( int hypId, SMESH_Gen* gen) + : StdMeshers_Regular_1D( hypId, gen) { } // ----------------------------------------------------------------------------- @@ -745,18 +774,30 @@ protected: * \brief Allow algo to do something after persistent restoration * \param subMesh - restored submesh * - * call markEdgeAsComputedByMe() + * call TEdgeMarker::markEdge() */ //======================================================================= void StdMeshers_RadialQuadrangle_1D2D::SubmeshRestored(SMESH_subMesh* faceSubMesh) { if ( !faceSubMesh->IsEmpty() ) + SetEventListener( faceSubMesh ); +} + +//======================================================================= +/*! + * \brief Sets event listener to a submesh + * \param subMesh - submesh where algo is set + * + * This method is called when a submesh gets HYP_OK algo_state. + */ +//======================================================================= + +void StdMeshers_RadialQuadrangle_1D2D::SetEventListener(SMESH_subMesh* faceSubMesh) +{ + for ( TopExp_Explorer e( faceSubMesh->GetSubShape(), TopAbs_EDGE ); e.More(); e.Next() ) { - for ( TopExp_Explorer e( faceSubMesh->GetSubShape(), TopAbs_EDGE ); e.More(); e.Next() ) - { - markEdgeAsComputedByMe( TopoDS::Edge( e.Current() ), faceSubMesh ); - } + TEdgeMarker::markEdge( TopoDS::Edge( e.Current() ), faceSubMesh ); } } @@ -964,7 +1005,7 @@ bool StdMeshers_RadialQuadrangle_1D2D::Compute(SMESH_Mesh& aMesh, list< TopoDS_Edge >::iterator ee = emptyEdges.begin(); for ( ; ee != emptyEdges.end(); ++ee ) - markEdgeAsComputedByMe( *ee, aMesh.GetSubMesh( F )); + TEdgeMarker::markEdge( *ee, aMesh.GetSubMesh( F )); circSide->GetUVPtStruct(); // let sides take into account just computed nodes linSide1->GetUVPtStruct(); @@ -976,6 +1017,7 @@ bool StdMeshers_RadialQuadrangle_1D2D::Compute(SMESH_Mesh& aMesh, quad->side[1] = linSide1; quad->side[2] = StdMeshers_FaceSide::New( circSide.get(), centerNode, ¢erUV ); quad->side[3] = linSide2; + quad->face = F; myQuadList.push_back( quad ); @@ -986,6 +1028,12 @@ bool StdMeshers_RadialQuadrangle_1D2D::Compute(SMESH_Mesh& aMesh, else ok = StdMeshers_Quadrangle_2D::computeTriangles( aMesh, F, quad ); + if ( helper.HasDegeneratedEdges() ) + { + StdMeshers_Quadrangle_2D::myNeedSmooth = true; + StdMeshers_Quadrangle_2D::smooth( quad ); + } + StdMeshers_Quadrangle_2D::myHelper = 0; return ok; @@ -1098,8 +1146,8 @@ bool StdMeshers_RadialQuadrangle_1D2D::Evaluate(SMESH_Mesh& aMesh, if( aResMap.count(sm) ) return false; - vector& aResVec = - aResMap.insert( make_pair(sm, vector(SMDSEntity_Last,0))).first->second; + vector& aResVec = + aResMap.insert( make_pair(sm, vector(SMDSEntity_Last,0))).first->second; myHelper = new SMESH_MesherHelper( aMesh ); myHelper->SetSubShape( aShape ); @@ -1156,23 +1204,23 @@ bool StdMeshers_RadialQuadrangle_1D2D::Evaluate(SMESH_Mesh& aMesh, for ( TopExp_Explorer edge( aShape, TopAbs_EDGE ); edge.More() && !isQuadratic ; edge.Next() ) { sm = aMesh.GetSubMesh( edge.Current() ); - vector& nbElems = aResMap[ sm ]; + vector& nbElems = aResMap[ sm ]; if ( SMDSEntity_Quad_Edge < (int) nbElems.size() ) isQuadratic = nbElems[ SMDSEntity_Quad_Edge ]; } - int nbCircSegments = 0; + smIdType nbCircSegments = 0; for ( int iE = 0; iE < circSide->NbEdges(); ++iE ) { sm = aMesh.GetSubMesh( circSide->Edge( iE )); - vector& nbElems = aResMap[ sm ]; + vector& nbElems = aResMap[ sm ]; if ( SMDSEntity_Quad_Edge < (int) nbElems.size() ) nbCircSegments += ( nbElems[ SMDSEntity_Edge ] + nbElems[ SMDSEntity_Quad_Edge ]); } - int nbQuads = nbCircSegments * ( layerPositions.size() - 1 ); - int nbTria = nbCircSegments; - int nbNodes = ( nbCircSegments - 1 ) * ( layerPositions.size() - 2 ); + smIdType nbQuads = nbCircSegments * ( layerPositions.size() - 1 ); + smIdType nbTria = nbCircSegments; + smIdType nbNodes = ( nbCircSegments - 1 ) * ( layerPositions.size() - 2 ); if ( isQuadratic ) { nbNodes += (( nbCircSegments - 1 ) * ( layerPositions.size() - 1 ) + // radial @@ -1190,7 +1238,7 @@ bool StdMeshers_RadialQuadrangle_1D2D::Evaluate(SMESH_Mesh& aMesh, if ( linSide1 ) { // evaluation for linSides - vector aResVec(SMDSEntity_Last, 0); + vector aResVec(SMDSEntity_Last, 0); if ( isQuadratic ) { aResVec[SMDSEntity_Node ] = 2 * ( layerPositions.size() - 1 ) + 1; aResVec[SMDSEntity_Quad_Edge] = layerPositions.size() - 1; @@ -1231,4 +1279,4 @@ bool StdMeshers_RadialQuadrangle_1D2D::IsApplicable( const TopoDS_Shape & aShape } if( toCheckAll && nbFoundFaces != 0 ) return true; return false; -}; +}