IPAL0052448: Too thin viscous layers are constructed

[modules/smesh.git] / src / StdMeshers / StdMeshers_ViscousLayers.cxx

diff --git a/src/StdMeshers/StdMeshers_ViscousLayers.cxx b/src/StdMeshers/StdMeshers_ViscousLayers.cxx
index 685b141fe09831619318af086650c4caebf9b845..241a21636aa89f76dc89d610accf9120ac302e8c 100644 (file)
--- a/src/StdMeshers/StdMeshers_ViscousLayers.cxx
+++ b/src/StdMeshers/StdMeshers_ViscousLayers.cxx
@@ -897,12 +897,12 @@ namespace
         if ( SMESH_Algo::isDegenerated( e )) continue;
         TopExp::Vertices( e, VV[0], VV[1], /*CumOri=*/true );
         if ( VV[1].IsSame( fromV )) {
+          nbEdges += edges[ 0 ].IsNull();
           edges[ 0 ] = e;
-          nbEdges++;
         }
         else if ( VV[0].IsSame( fromV )) {
+          nbEdges += edges[ 1 ].IsNull();
           edges[ 1 ] = e;
-          nbEdges++;
         }
       }
     }
@@ -922,7 +922,7 @@ namespace
 
       // get angle between the 2 edges
       gp_Vec faceNormal;
-      double angle = helper.GetAngle( edges[0], edges[1], faceFrw, &faceNormal );
+      double angle = helper.GetAngle( edges[0], edges[1], faceFrw, fromV, &faceNormal );
       if ( Abs( angle ) < 5 * M_PI/180 )
       {
         dir = ( faceNormal.XYZ() ^ edgeDir[0].Reversed()) + ( faceNormal.XYZ() ^ edgeDir[1] );
@@ -940,7 +940,7 @@ namespace
     }
     else if ( nbEdges == 1 )
     {
-      dir = getFaceDir( faceFrw, edges[0], node, helper, ok );
+      dir = getFaceDir( faceFrw, edges[ edges[0].IsNull() ], node, helper, ok );
       if ( cosin ) *cosin = 1.;
     }
     else
@@ -1004,7 +1004,8 @@ namespace
         while ( SMESH_Algo::isDegenerated( wires[iW]->Edge( iE2 )))
           iE2 = ( iE2 + 1 ) % nbEdges;
         double angle = helper.GetAngle( wires[iW]->Edge( iE1 ),
-                                        wires[iW]->Edge( iE2 ), F );
+                                        wires[iW]->Edge( iE2 ), F,
+                                        wires[iW]->FirstVertex( iE2 ));
         if ( angle < -5. * M_PI / 180. )
           return true;
       }
@@ -2283,12 +2284,19 @@ gp_XYZ _ViscousBuilder::getFaceNormal(const SMDS_MeshNode* node,
   isOK = false;
 
   Handle(Geom_Surface) surface = BRep_Tool::Surface( face );
-  if ( GeomLib::NormEstim( surface, uv, 1e-10, normal ) < 3 )
+  int pointKind = GeomLib::NormEstim( surface, uv, 1e-5, normal );
+  enum { REGULAR = 0, QUASYSINGULAR, CONICAL, IMPOSSIBLE };
+  if ( pointKind < IMPOSSIBLE )
   {
-    normal;
+    if ( pointKind != REGULAR && !shiftInside )
+    {
+      gp_XYZ normShift = getFaceNormal( node, face, helper, isOK, /*shiftInside=*/true );
+      if ( normShift * normal.XYZ() < 0. )
+        normal = normShift;
+    }
     isOK = true;
   }
-  else // hard singularity
+  else // hard singularity, to call with shiftInside=true ?
   {
     const TGeomID faceID = helper.GetMeshDS()->ShapeToIndex( face );
 
@@ -2371,12 +2379,10 @@ gp_XYZ _ViscousBuilder::getWeigthedNormal( const SMDS_MeshNode*         n,
     }
     else
     {
-      TopoDS_Vertex v10 = SMESH_MesherHelper::IthVertex( 1, ee[ 0 ]);
-      TopoDS_Vertex v01 = SMESH_MesherHelper::IthVertex( 0, ee[ 1 ]);
-      if ( !v10.IsSame( v01 ))
+      if ( !V.IsSame( SMESH_MesherHelper::IthVertex( 0, ee[ 1 ] )))
         std::swap( ee[0], ee[1] );
     }
-    angles[i] = SMESH_MesherHelper::GetAngle( ee[0], ee[1], F );
+    angles[i] = SMESH_MesherHelper::GetAngle( ee[0], ee[1], F, TopoDS::Vertex( V ));
   }
 
   // compute a weighted normal