-// Copyright (C) 2007-2021 CEA/DEN, EDF R&D, OPEN CASCADE
+// Copyright (C) 2007-2024 CEA, EDF, OPEN CASCADE
//
// This library is free software; you can redistribute it and/or
// modify it under the terms of the GNU Lesser General Public
#include "StdMeshers_Quadrangle_2D.hxx"
#include "StdMeshers_ViscousLayers2D.hxx"
+#include <Basics_OCCTVersion.hxx>
+
+#if OCC_VERSION_LARGE < 0x07070000
#include <Adaptor3d_HSurface.hxx>
+#else
+#include <Adaptor3d_Surface.hxx>
+#endif
+
#include <BRepAdaptor_Curve.hxx>
#include <BRepAdaptor_Curve2d.hxx>
#include <BRepAdaptor_Surface.hxx>
#include <unordered_map>
#ifdef _DEBUG_
+#ifndef WIN32
#define __myDEBUG
+#endif
//#define __NOT_INVALIDATE_BAD_SMOOTH
//#define __NODES_AT_POS
#endif
const double T = ( realThickness > 0 ) ? realThickness : GetTotalThickness();
const double f = GetStretchFactor();
const int N = GetNumberLayers();
- const double fPowN = pow( f, N );
- double h0;
- if ( fPowN - 1 <= numeric_limits<double>::min() )
- h0 = T / N;
- else
- h0 = T * ( f - 1 )/( fPowN - 1 );
- return h0;
+ return StdMeshers_ViscousLayers::Get1stLayerThickness( T, f, N );
}
bool UseSurfaceNormal() const
( std::find( _shapeIds.begin(), _shapeIds.end(), shapeIndex ) != _shapeIds.end() );
return IsToIgnoreShapes() ? !isIn : isIn;
}
-
+// --------------------------------------------------------------------------------
+double StdMeshers_ViscousLayers::Get1stLayerThickness( double T, double f, int N )
+{
+ const double fPowN = pow( f, N );
+ double h0;
+ if ( fPowN - 1 <= numeric_limits<double>::min() )
+ h0 = T / N;
+ else
+ h0 = T * ( f - 1 )/( fPowN - 1 );
+ return h0;
+}
// --------------------------------------------------------------------------------
SMDS_MeshGroup* StdMeshers_ViscousLayers::CreateGroup( const std::string& theName,
SMESH_Mesh& theMesh,
//case GeomAbs_SurfaceOfExtrusion:
case GeomAbs_OffsetSurface:
{
+#if OCC_VERSION_LARGE < 0x07070000
Handle(Adaptor3d_HSurface) base = surface.BasisSurface();
return getRovolutionAxis( base->Surface(), axis );
+#else
+ Handle(Adaptor3d_Surface) base = surface.BasisSurface();
+ return getRovolutionAxis( *base, axis );
+#endif
}
default: return false;
}
if ( ! shrink(_sdVec[iSD]) ) // shrink 2D mesh on FACEs w/o layer
return _error;
- addBoundaryElements(_sdVec[iSD]); // create quadrangles on prism bare sides
+ bool notMissingFaces = addBoundaryElements(_sdVec[iSD]); // create quadrangles on prism bare sides
+ if ( !notMissingFaces )
+ {
+ SMESH_MeshEditor editor( &theMesh );
+ TIDSortedElemSet elements;
+ editor.MakeBoundaryMesh( elements, SMESH_MeshEditor::BND_2DFROM3D );
+ }
+
_sdVec[iSD]._done = true;
const TopoDS_Shape& solid = _sdVec[iSD]._solid;
break;
}
default:
- return error("Not yet supported case", _sdVec[i]._index);
+ std::ostringstream msg;
+ msg << "Not yet supported case: vertex bounded by ";
+ msg << facesWOL.size();
+ msg << " faces without layer at coordinates (";
+ TopoDS_Vertex v = TopoDS::Vertex(vertex);
+ gp_Pnt p = BRep_Tool::Pnt(v);
+ msg << p.X() << ", " << p.Y() << ", " << p.Z() << ")";
+ return error(msg.str().c_str(), _sdVec[i]._index);
}
}
}
if ( !setEdgeData( *edge, edgesByGeom[ shapeID ], helper, data ))
return false;
- if ( edge->_nodes.size() < 2 )
- edge->Block( data );
- //data._noShrinkShapes.insert( shapeID );
+ if ( edge->_nodes.size() < 2 && !noShrink )
+ edge->Block( data ); // a sole node is moved only if noShrink
}
dumpMove(edge->_nodes.back());
- if ( edge->_cosin > faceMaxCosin && !edge->Is( _LayerEdge::BLOCKED ))
+ if ( edge->_cosin > faceMaxCosin && edge->_nodes.size() > 1 )
{
faceMaxCosin = edge->_cosin;
maxCosinEdge = edge;
cnvFace._face = F;
cnvFace._normalsFixed = false;
cnvFace._isTooCurved = false;
+ cnvFace._normalsFixedOnBorders = false;
double maxCurvature = cnvFace.GetMaxCurvature( data, eof, surfProp, helper );
if ( maxCurvature > 0 )
// Find C1 EDGEs
vector< pair< _EdgesOnShape*, gp_XYZ > > dirOfEdges;
-
+
for ( size_t iS = 0; iS < edgesByGeom.size(); ++iS ) // check VERTEXes
{
_EdgesOnShape& eov = edgesByGeom[iS];
if ( eov._edges.empty() ||
eov.ShapeType() != TopAbs_VERTEX ||
c1VV.Contains( eov._shape ))
- continue;
+ continue;
const TopoDS_Vertex& V = TopoDS::Vertex( eov._shape );
// get directions of surrounding EDGEs
if ( oppV.IsSame( V ))
oppV = SMESH_MesherHelper::IthVertex( 1, e );
_EdgesOnShape* eovOpp = data.GetShapeEdges( oppV );
- if ( dirOfEdges[k].second * eovOpp->_edges[0]->_normal < 0 )
+ if ( !eovOpp->_edges.empty() && dirOfEdges[k].second * eovOpp->_edges[0]->_normal < 0 )
eov._eosC1.push_back( dirOfEdges[k].first );
}
dirOfEdges[k].first = 0;
getMeshDS()->SetNodeOnFace( tgtNode, TopoDS::Face( eos._sWOL ), uv.X(), uv.Y() );
}
- if ( edge._nodes.size() > 1 )
+ //if ( edge._nodes.size() > 1 ) -- allow RISKY_SWOL on noShrink shape
{
// check if an angle between a FACE with layers and SWOL is sharp,
// else the edge should not inflate
geomNorm.Reverse(); // inside the SOLID
if ( geomNorm * edge._normal < -0.001 )
{
- getMeshDS()->RemoveFreeNode( tgtNode, 0, /*fromGroups=*/false );
- edge._nodes.resize( 1 );
+ if ( edge._nodes.size() > 1 )
+ {
+ getMeshDS()->RemoveFreeNode( tgtNode, 0, /*fromGroups=*/false );
+ edge._nodes.resize( 1 );
+ }
}
else if ( realLenFactor > 3 ) /// -- moved to SetCosin()
//else if ( edge._lenFactor > 3 )
void _ViscousBuilder::makeGroupOfLE()
{
-#ifdef _DEBUG_
+ if (!SALOME::VerbosityActivated())
+ return;
+
for ( size_t i = 0 ; i < _sdVec.size(); ++i )
{
if ( _sdVec[i]._n2eMap.empty() ) continue;
<< "'%s-%s' % (faceId1+1, faceId2))");
dumpFunctionEnd();
}
-#endif
}
//================================================================================
double thinkness = eos._hyp.GetTotalThickness();
for ( size_t i = 0; i < eos._edges.size(); ++i )
{
- if ( eos._edges[i]->Is( _LayerEdge::BLOCKED )) continue;
+ if ( eos._edges[i]->_nodes.size() < 2 ) continue;
eos._edges[i]->SetMaxLen( thinkness );
eos._edges[i]->FindIntersection( *searcher, intersecDist, data._epsilon, eos, &face );
if ( intersecDist > 0 && face )
if ( eos._edges[i]->_nodes.size() > 1 )
avgThick += Min( 1., eos._edges[i]->_len / shapeTgtThick );
else
- avgThick += shapeTgtThick;
+ avgThick += 1;
nbActiveEdges += ( ! eos._edges[i]->Is( _LayerEdge::BLOCKED ));
}
}
// intersection not ignored
- if ( toBlockInfaltion &&
- dist < ( eos._edges[i]->_len * theThickToIntersection ))
+ double minDist = 0;
+ if ( eos._edges[i]->_maxLen < 0.99 * eos._hyp.GetTotalThickness() ) // limited length
+ minDist = eos._edges[i]->_len * theThickToIntersection;
+
+ if ( toBlockInfaltion && dist < minDist )
{
if ( is1stBlocked ) { is1stBlocked = false; // debug
dumpFunction(SMESH_Comment("blockIntersected") <<data._index<<"_InfStep"<<infStep);
}
}
-
-
-#ifdef _DEBUG_
- // dumpMove() for debug
- size_t i = 0;
- for ( ; i < eos._edges.size(); ++i )
- if ( eos._edges[i]->Is( _LayerEdge::MARKED ))
- break;
- if ( i < eos._edges.size() )
+ if (SALOME::VerbosityActivated())
{
- dumpFunction(SMESH_Comment("putOnOffsetSurface_") << eos.ShapeTypeLetter() << eos._shapeID
- << "_InfStep" << infStep << "_" << Abs( smooStep ));
+ // dumpMove() for debug
+ size_t i = 0;
for ( ; i < eos._edges.size(); ++i )
+ if ( eos._edges[i]->Is( _LayerEdge::MARKED ))
+ break;
+ if ( i < eos._edges.size() )
{
- if ( eos._edges[i]->Is( _LayerEdge::MARKED )) {
- dumpMove( eos._edges[i]->_nodes.back() );
+ dumpFunction(SMESH_Comment("putOnOffsetSurface_") << eos.ShapeTypeLetter() << eos._shapeID
+ << "_InfStep" << infStep << "_" << Abs( smooStep ));
+ for ( ; i < eos._edges.size(); ++i )
+ {
+ if ( eos._edges[i]->Is( _LayerEdge::MARKED )) {
+ dumpMove( eos._edges[i]->_nodes.back() );
+ }
}
+ dumpFunctionEnd();
}
- dumpFunctionEnd();
}
-#endif
_ConvexFace* cnvFace;
if ( moveAll != _LayerEdge::UPD_NORMAL_CONV &&
tgtNode->setXYZ( newPos.X(), newPos.Y(), newPos.Z() );
dumpMove( tgtNode );
- SMDS_FacePositionPtr pos = tgtNode->GetPosition();
- pos->SetUParameter( newUV.X() );
- pos->SetVParameter( newUV.Y() );
+ if ( SMDS_FacePositionPtr pos = tgtNode->GetPosition() ) // NULL if F is noShrink
+ {
+ pos->SetUParameter( newUV.X() );
+ pos->SetVParameter( newUV.Y() );
+ }
gp_XYZ newUV0( newUV.X(), newUV.Y(), 0 );
tgtNode->setXYZ( newPos.X(), newPos.Y(), newPos.Z() );
dumpMove( tgtNode );
- SMDS_FacePositionPtr pos = tgtNode->GetPosition();
- pos->SetUParameter( newUV.X() );
- pos->SetVParameter( newUV.Y() );
-
+ if ( SMDS_FacePositionPtr pos = tgtNode->GetPosition() ) // NULL if F is noShrink
+ {
+ pos->SetUParameter( newUV.X() );
+ pos->SetVParameter( newUV.Y() );
+ }
_eos[i]->Set( _LayerEdge::SMOOTHED ); // to check in refine() (IPAL54237)
}
}
_edgeDir[1] = getEdgeDir( E, leOnV[1]->_nodes[0], data.GetHelper() );
_leOnV[0]._cosin = Abs( leOnV[0]->_cosin );
_leOnV[1]._cosin = Abs( leOnV[1]->_cosin );
+ _leOnV[0]._flags = _leOnV[1]._flags = 0;
if ( _eos._sWOL.IsNull() ) // 3D
for ( int iEnd = 0; iEnd < 2; ++iEnd )
_leOnV[iEnd]._cosin = Abs( _edgeDir[iEnd].Normalized() * leOnV[iEnd]->_normal );
eos->_edgeForOffset = 0;
double maxCosin = -1;
- bool hasNoShrink = false;
+ //bool hasNoShrink = false;
for ( TopExp_Explorer eExp( eos->_shape, TopAbs_EDGE ); eExp.More(); eExp.Next() )
{
_EdgesOnShape* eoe = GetShapeEdges( eExp.Current() );
if ( !eoe || eoe->_edges.empty() ) continue;
- if ( eos->GetData()._noShrinkShapes.count( eoe->_shapeID ))
- hasNoShrink = true;
+ // if ( eos->GetData()._noShrinkShapes.count( eoe->_shapeID ))
+ // hasNoShrink = true;
vector<_LayerEdge*>& eE = eoe->_edges;
_LayerEdge* e = eE[ eE.size() / 2 ];
// Try to initialize _Mapper2D
- if ( hasNoShrink )
- return;
+ // if ( hasNoShrink )
+ // return;
SMDS_ElemIteratorPtr fIt = eos->_subMesh->GetSubMeshDS()->GetElements();
if ( !fIt->more() || fIt->next()->NbCornerNodes() != 4 )
std::vector< SMESH_NodeXYZ > _nodes;
TopAbs_ShapeEnum _vertSWOLType[2]; // shrink part includes VERTEXes
AverageHyp* _vertHyp[2];
+ double _edgeWOLLen[2]; // length of wol EDGE
+ double _tol; // to compare _edgeWOLLen's
BndPart():
_isShrink(0), _isReverse(0), _nbSegments(0), _hyp(0),
- _vertSWOLType{ TopAbs_WIRE, TopAbs_WIRE }, _vertHyp{ 0, 0 }
+ _vertSWOLType{ TopAbs_WIRE, TopAbs_WIRE }, _vertHyp{ 0, 0 }, _edgeWOLLen{ 0., 0.}
{}
+ bool IsEqualLengthEWOL( const BndPart& other ) const
+ {
+ return ( std::abs( _edgeWOLLen[0] - other._edgeWOLLen[0] ) < _tol &&
+ std::abs( _edgeWOLLen[1] - other._edgeWOLLen[1] ) < _tol );
+ }
+
bool operator==( const BndPart& other ) const
{
return ( _isShrink == other._isShrink &&
(( !_isShrink ) ||
( *_hyp == *other._hyp &&
vertHyp1() == other.vertHyp1() &&
- vertHyp2() == other.vertHyp2() ))
+ vertHyp2() == other.vertHyp2() &&
+ IsEqualLengthEWOL( other )))
);
}
bool CanAppend( const BndPart& other )
bool hasCommonNode = ( _nodes.back()->GetID() == other._nodes.front()->GetID() );
_nodes.insert( _nodes.end(), other._nodes.begin() + hasCommonNode, other._nodes.end() );
_vertSWOLType[1] = other._vertSWOLType[1];
- if ( _isShrink )
- _vertHyp[1] = other._vertHyp[1];
+ if ( _isShrink ) {
+ _vertHyp[1] = other._vertHyp[1];
+ _edgeWOLLen[1] = other._edgeWOLLen[1];
+ }
}
- const SMDS_MeshNode* Node(size_t i) const
+ const SMDS_MeshNode* Node(size_t i) const
{
return _nodes[ _isReverse ? ( _nodes.size() - 1 - i ) : i ]._node;
}
for ( int iE = 0; iE < nbEdgesInWire.front(); ++iE )
{
BndPart bndPart;
+
+ std::vector<const SMDS_MeshNode*> nodes = fSide.GetOrderedNodes( iE );
+ bndPart._nodes.assign( nodes.begin(), nodes.end() );
+ bndPart._nbSegments = bndPart._nodes.size() - 1;
+
_EdgesOnShape* eos = _data1->GetShapeEdges( fSide.EdgeID( iE ));
bndPart._isShrink = ( eos->SWOLType() == TopAbs_FACE );
bndPart._vertSWOLType[iV] = eov[iV]->SWOLType();
}
}
+ bndPart._edgeWOLLen[0] = fSide.EdgeLength( iE - 1 );
+ bndPart._edgeWOLLen[1] = fSide.EdgeLength( iE + 1 );
+
+ bndPart._tol = std::numeric_limits<double>::max(); // tolerance by segment size
+ for ( size_t i = 1; i < bndPart._nodes.size(); ++i )
+ bndPart._tol = Min( bndPart._tol,
+ ( bndPart._nodes[i-1] - bndPart._nodes[i] ).SquareModulus() );
}
- std::vector<const SMDS_MeshNode*> nodes = fSide.GetOrderedNodes( iE );
- bndPart._nodes.assign( nodes.begin(), nodes.end() );
- bndPart._nbSegments = bndPart._nodes.size() - 1;
if ( _boundary.empty() || ! _boundary.back().CanAppend( bndPart ))
_boundary.push_back( bndPart );
}
SMESHDS_Mesh* meshDS = dataSrc->GetHelper().GetMeshDS();
+ dumpFunction(SMESH_Comment("periodicMoveNodes_F")
+ << _shriFace[iSrc]->_subMesh->GetId() << "_F"
+ << _shriFace[iTgt]->_subMesh->GetId() );
TNode2Edge::iterator n2e;
TNodeNodeMap::iterator n2n = _nnMap.begin();
for ( ; n2n != _nnMap.end(); ++n2n )
SMESH_NodeXYZ pSrc = leSrc->_nodes[ iN ];
gp_XYZ pTgt = trsf->Transform( pSrc );
meshDS->MoveNode( leTgt->_nodes[ iN ], pTgt.X(), pTgt.Y(), pTgt.Z() );
+
+ dumpMove( leTgt->_nodes[ iN ]);
}
}
}
<< _shriFace[iSrc]->_subMesh->GetId() << " -> "
<< _shriFace[iTgt]->_subMesh->GetId() << " -- "
<< ( done ? "DONE" : "FAIL"));
+ dumpFunctionEnd();
return done;
}
{
_EdgesOnShape& eos = * subEOS[ iS ];
if ( eos.ShapeType() != TopAbs_EDGE ) continue;
+ if ( eos.size() == 0 )
+ continue;
const TopoDS_Edge& E = TopoDS::Edge( eos._shape );
data.SortOnEdge( E, eos._edges );
uvPtVec[ i ].param = helper.GetNodeU( E, edges[i]->_nodes[0] );
uvPtVec[ i ].SetUV( helper.GetNodeUV( F, edges[i]->_nodes.back() ));
}
- if ( edges.empty() )
- continue;
+ if ( uvPtVec[ 0 ].node == uvPtVec.back().node && // closed
+ helper.IsSeamShape( uvPtVec[ 0 ].node->GetShapeID() ))
+ {
+ uvPtVec[ 0 ].SetUV( helper.GetNodeUV( F,
+ edges[0]->_nodes.back(),
+ edges[1]->_nodes.back() ));
+ size_t i = edges.size() - 1;
+ uvPtVec[ i ].SetUV( helper.GetNodeUV( F,
+ edges[i ]->_nodes.back(),
+ edges[i-1]->_nodes.back() ));
+ }
+ // if ( edges.empty() )
+ // continue;
BRep_Tool::Range( E, uvPtVec[0].param, uvPtVec.back().param );
StdMeshers_FaceSide fSide( uvPtVec, F, E, _mesh );
StdMeshers_ViscousLayers2D::SetProxyMeshOfEdge( fSide );
smDS->Clear();
// compute the mesh on the FACE
+ TopTools_IndexedMapOfShape allowed(1);
+ allowed.Add( sm->GetSubShape() );
+ sm->SetAllowedSubShapes( &allowed );
sm->ComputeStateEngine( SMESH_subMesh::CHECK_COMPUTE_STATE );
sm->ComputeStateEngine( SMESH_subMesh::COMPUTE_SUBMESH );
+ sm->SetAllowedSubShapes( nullptr );
// re-fill proxy sub-meshes of the FACE
for ( size_t i = 0 ; i < _sdVec.size(); ++i )
edgeSize.back() = edgeSize.front();
gp_XY newPos(0,0);
- //int nbEdges = 0;
- double sumSize = 0;
+ double sumWgt = 0;
for ( size_t i = 1; i < edgeDir.size(); ++i )
{
- if ( edgeDir[i-1].X() > 1. ) continue;
- int i1 = i-1;
+ const int i1 = i-1;
+ if ( edgeDir[i1].X() > 1. ) continue;
while ( edgeDir[i].X() > 1. && ++i < edgeDir.size() );
if ( i == edgeDir.size() ) break;
gp_XY p = uv[i];
gp_XY norm1( -edgeDir[i1].Y(), edgeDir[i1].X() );
- gp_XY norm2( -edgeDir[i].Y(), edgeDir[i].X() );
+ gp_XY norm2( -edgeDir[i ].Y(), edgeDir[i ].X() );
gp_XY bisec = norm1 + norm2;
double bisecSize = bisec.Modulus();
if ( bisecSize < numeric_limits<double>::min() )
}
bisec /= bisecSize;
- gp_XY dirToN = uvToFix - p;
- double distToN = dirToN.Modulus();
+ gp_XY dirToN = uvToFix - p;
+ double distToN = bisec * dirToN;
if ( bisec * dirToN < 0 )
distToN = -distToN;
- newPos += ( p + bisec * distToN ) * ( edgeSize[i1] + edgeSize[i] );
- //++nbEdges;
- sumSize += edgeSize[i1] + edgeSize[i];
+ double wgt = edgeSize[i1] + edgeSize[i];
+ newPos += ( p + bisec * distToN ) * wgt;
+ sumWgt += wgt;
}
- newPos /= /*nbEdges * */sumSize;
+ newPos /= sumWgt;
return newPos;
}
bool _ViscousBuilder::addBoundaryElements(_SolidData& data)
{
+ bool addAllBoundaryElements = true;
SMESH_MesherHelper helper( *_mesh );
vector< const SMDS_MeshNode* > faceNodes;
map< double, const SMDS_MeshNode* > u2nodes;
if ( !SMESH_Algo::GetSortedNodesOnEdge( getMeshDS(), E, /*ignoreMedium=*/false, u2nodes))
continue;
-
+
vector< _LayerEdge* > ledges; ledges.reserve( u2nodes.size() );
TNode2Edge & n2eMap = data._n2eMap;
map< double, const SMDS_MeshNode* >::iterator u2n = u2nodes.begin();
if ( getMeshDS()->FindElement( faceNodes, SMDSAbs_Face, /*noMedium=*/true))
continue; // faces already created
- }
+ }
for ( ++u2n; u2n != u2nodes.end(); ++u2n )
- ledges.push_back( n2eMap[ u2n->second ]);
+ if ( n2eMap[ u2n->second ] != nullptr )
+ ledges.push_back( n2eMap[ u2n->second ]);
+ else /*some boundary elements might be lost because the connectivity of the face is not entirely defined on this edge*/
+ addAllBoundaryElements = false;
// Find out orientation and type of face to create
-
bool reverse = false, isOnFace;
TopoDS_Shape F;
} // loop on EDGE's
} // loop on _SolidData's
- return true;
+ return addAllBoundaryElements;
}
