Salome HOME
Cleanup of parallel meshing + documentation
[modules/smesh.git] / src / StdMeshers / StdMeshers_ViscousLayers.cxx
index ee90bca54f66b1e62276af0889e5cefd0c43129e..a8d421abaffe3623b76d2e6171d5168497995deb 100644 (file)
@@ -1,4 +1,4 @@
-// Copyright (C) 2007-2021  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
 #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
@@ -1816,8 +1825,13 @@ namespace VISCOUS_3D
     //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;
     }
@@ -2429,7 +2443,14 @@ bool _ViscousBuilder::findFacesWithLayers(const bool onlyWith)
         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);
       }
     }
   }
@@ -3036,6 +3057,7 @@ void _ViscousBuilder::limitStepSizeByCurvature( _SolidData& data )
     cnvFace._face = F;
     cnvFace._normalsFixed = false;
     cnvFace._isTooCurved = false;
+    cnvFace._normalsFixedOnBorders = false;
 
     double maxCurvature = cnvFace.GetMaxCurvature( data, eof, surfProp, helper );
     if ( maxCurvature > 0 )
@@ -4641,7 +4663,9 @@ void _Simplex::SortSimplices(vector<_Simplex>& simplices)
 
 void _ViscousBuilder::makeGroupOfLE()
 {
-#ifdef _DEBUG_
+  if (!SALOME::VerbosityActivated())
+    return;
+
   for ( size_t i = 0 ; i < _sdVec.size(); ++i )
   {
     if ( _sdVec[i]._n2eMap.empty() ) continue;
@@ -4697,7 +4721,6 @@ void _ViscousBuilder::makeGroupOfLE()
              << "'%s-%s' % (faceId1+1, faceId2))");
     dumpFunctionEnd();
   }
-#endif
 }
 
 //================================================================================
@@ -4932,7 +4955,7 @@ bool _ViscousBuilder::inflate(_SolidData& data)
         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 ));
       }
     }
@@ -5446,8 +5469,11 @@ bool _ViscousBuilder::smoothAndCheck(_SolidData& data,
 
         // 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);
@@ -5829,27 +5855,26 @@ void _ViscousBuilder::putOnOffsetSurface( _EdgesOnShape&            eos,
     }
   }
 
-
-
-#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 &&
@@ -6560,6 +6585,7 @@ void _Smoother1D::prepare(_SolidData& data)
   _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 );
@@ -11790,6 +11816,17 @@ bool _ViscousBuilder::shrink(_SolidData& theData)
           uvPtVec[ i ].param = helper.GetNodeU( E, edges[i]->_nodes[0] );
           uvPtVec[ i ].SetUV( helper.GetNodeUV( F, edges[i]->_nodes.back() ));
         }
+        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 );
@@ -11806,8 +11843,12 @@ bool _ViscousBuilder::shrink(_SolidData& theData)
       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 )
@@ -12408,17 +12449,16 @@ gp_XY _SmoothNode::computeAngularPos(vector<gp_XY>& uv,
   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() )
@@ -12428,16 +12468,16 @@ gp_XY _SmoothNode::computeAngularPos(vector<gp_XY>& uv,
     }
     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;
 }