Salome HOME
Fix regressions
[modules/smesh.git] / src / StdMeshers / StdMeshers_ViscousLayers.cxx
index f84fc5b5e791102a8d181afa598a336bcce75639..a6530b9244a5782abe3a820e82075e2aab951077 100644 (file)
@@ -95,7 +95,7 @@
 #include <string>
 
 #ifdef _DEBUG_
-#define __myDEBUG
+//#define __myDEBUG
 //#define __NOT_INVALIDATE_BAD_SMOOTH
 //#define __NODES_AT_POS
 #endif
@@ -409,7 +409,7 @@ namespace VISCOUS_3D
 
     gp_XYZ              _normal;    // to boundary of solid
     vector<gp_XYZ>      _pos;       // points computed during inflation
-    double              _len;       // length achived with the last inflation step
+    double              _len;       // length achieved with the last inflation step
     double              _maxLen;    // maximal possible length
     double              _cosin;     // of angle (_normal ^ surface)
     double              _minAngle;  // of _simplices
@@ -440,11 +440,12 @@ namespace VISCOUS_3D
                   DISTORTED       = 0x0001000, // was bad before smoothing
                   RISKY_SWOL      = 0x0002000, // SWOL is parallel to a source FACE
                   SHRUNK          = 0x0004000, // target node reached a tgt position while shrink()
-                  UNUSED_FLAG     = 0x0100000
+                  UNUSED_FLAG     = 0x0100000  // to add use flags after
     };
     bool Is   ( int flag ) const { return _flags & flag; }
     void Set  ( int flag ) { _flags |= flag; }
     void Unset( int flag ) { _flags &= ~flag; }
+    std::string DumpFlags() const; // debug
 
     void SetNewLength( double len, _EdgesOnShape& eos, SMESH_MesherHelper& helper );
     bool SetNewLength2d( Handle(Geom_Surface)& surface,
@@ -1864,15 +1865,19 @@ SMESH_ComputeErrorPtr _ViscousBuilder::Compute(SMESH_Mesh&         theMesh,
     size_t iSD = 0;
     for ( iSD = 0; iSD < _sdVec.size(); ++iSD ) // find next SOLID to compute
       if ( _sdVec[iSD]._before.IsEmpty() &&
+           !_sdVec[iSD]._solid.IsNull() &&
            _sdVec[iSD]._n2eMap.empty() )
         break;
 
     if ( ! makeLayer(_sdVec[iSD]) )   // create _LayerEdge's
       return _error;
 
-    if ( _sdVec[iSD]._n2eMap.size() == 0 )
+    if ( _sdVec[iSD]._n2eMap.size() == 0 ) // no layers in a SOLID
+    {
+      _sdVec[iSD]._solid.Nullify();
       continue;
-    
+    }
+
     if ( ! inflate(_sdVec[iSD]) )     // increase length of _LayerEdge's
       return _error;
 
@@ -2948,7 +2953,7 @@ bool _ViscousBuilder::findShapesToSmooth( _SolidData& data )
   //const double tgtThick = /*Min( 0.5 * data._geomSize, */data._maxThickness;
 
   // Find shapes needing smoothing; such a shape has _LayerEdge._normal on it's
-  // boundry inclined to the shape at a sharp angle
+  // boundary inclined to the shape at a sharp angle
 
   //list< TGeomID > shapesToSmooth;
   TopTools_MapOfShape edgesOfSmooFaces;
@@ -4157,7 +4162,7 @@ void _LayerEdge::SetDataByNeighbors( const SMDS_MeshNode* n1,
 //================================================================================
 /*!
  * \brief Copy data from a _LayerEdge of other SOLID and based on the same node;
- * this and other _LayerEdge's are inflated along a FACE or an EDGE
+ * this and the other _LayerEdge are inflated along a FACE or an EDGE
  */
 //================================================================================
 
@@ -4507,7 +4512,10 @@ bool _ViscousBuilder::inflate(_SolidData& data)
       const double shapeTgtThick = eos._hyp.GetTotalThickness();
       for ( size_t i = 0; i < eos._edges.size(); ++i )
       {
-        avgThick      += Min( 1., eos._edges[i]->_len / shapeTgtThick );
+        if ( eos._edges[i]->_nodes.size() > 1 )
+          avgThick    += Min( 1., eos._edges[i]->_len / shapeTgtThick );
+        else
+          avgThick    += shapeTgtThick;
         nbActiveEdges += ( ! eos._edges[i]->Is( _LayerEdge::BLOCKED ));
       }
     }
@@ -8507,7 +8515,7 @@ gp_XYZ _LayerEdge::smoothNefPolygon()
 {
   gp_XYZ newPos(0,0,0);
 
-  // get a plane to seach a solution on
+  // get a plane to search a solution on
 
   vector< gp_XYZ > vecs( _simplices.size() + 1 );
   size_t i;
@@ -8698,7 +8706,7 @@ gp_XYZ _LayerEdge::smoothNefPolygon()
 { ////////////////////////////////// NEW
   gp_XYZ newPos(0,0,0);
 
-  // get a plane to seach a solution on
+  // get a plane to search a solution on
 
   size_t i;
   gp_XYZ center(0,0,0);
@@ -9045,7 +9053,7 @@ void _LayerEdge::Block( _SolidData& data )
       minDist   = Min( pSrc.SquareDistance( pTgtN ), minDist );
       minDist   = Min( pTgt.SquareDistance( pSrcN ), minDist );
       double newMaxLen = edge->_maxLen + 0.5 * Sqrt( minDist );
-      if ( edge->_nodes[0]->getshapeId() == neibor->_nodes[0]->getshapeId() )
+      //if ( edge->_nodes[0]->getshapeId() == neibor->_nodes[0]->getshapeId() ) viscous_layers_00/A3
       {
         newMaxLen *= edge->_lenFactor / neibor->_lenFactor;
       }
@@ -9199,10 +9207,48 @@ void _LayerEdge::SmoothPos( const vector< double >& segLen, const double tol )
 
 //================================================================================
 /*!
- * \brief Create layers of prisms
+ * \brief Print flags
  */
 //================================================================================
 
+std::string _LayerEdge::DumpFlags() const
+{
+  SMESH_Comment dump;
+  for ( int flag = 1; flag < 0x1000000; flag *= 2 )
+    if ( _flags & flag )
+    {
+      EFlags f = (EFlags) flag;
+      switch ( f ) {
+      case TO_SMOOTH:       dump << "TO_SMOOTH";       break;
+      case MOVED:           dump << "MOVED";           break;
+      case SMOOTHED:        dump << "SMOOTHED";        break;
+      case DIFFICULT:       dump << "DIFFICULT";       break;
+      case ON_CONCAVE_FACE: dump << "ON_CONCAVE_FACE"; break;
+      case BLOCKED:         dump << "BLOCKED";         break;
+      case INTERSECTED:     dump << "INTERSECTED";     break;
+      case NORMAL_UPDATED:  dump << "NORMAL_UPDATED";  break;
+      case MARKED:          dump << "MARKED";          break;
+      case MULTI_NORMAL:    dump << "MULTI_NORMAL";    break;
+      case NEAR_BOUNDARY:   dump << "NEAR_BOUNDARY";   break;
+      case SMOOTHED_C1:     dump << "SMOOTHED_C1";     break;
+      case DISTORTED:       dump << "DISTORTED";       break;
+      case RISKY_SWOL:      dump << "RISKY_SWOL";      break;
+      case SHRUNK:          dump << "SHRUNK";          break;
+      case UNUSED_FLAG:     dump << "UNUSED_FLAG";     break;
+      }
+      dump << " ";
+    }
+  cout << dump << endl;
+  return dump;
+}
+
+//================================================================================
+/*!
+  case brief:
+  default:
+*/
+//================================================================================
+
 bool _ViscousBuilder::refine(_SolidData& data)
 {
   SMESH_MesherHelper& helper = data.GetHelper();
@@ -9765,6 +9811,7 @@ bool _ViscousBuilder::shrink(_SolidData& theData)
       continue;
 
     _SolidData&      data = *dataList.front();
+    _SolidData*     data2 = dataList.size() > 1 ? dataList.back() : 0;
     const TopoDS_Face&  F = TopoDS::Face( getMeshDS()->IndexToShape( f2sd->first ));
     SMESH_subMesh*     sm = _mesh->GetSubMesh( F );
     SMESHDS_SubMesh* smDS = sm->GetSubMeshDS();
@@ -9778,7 +9825,7 @@ bool _ViscousBuilder::shrink(_SolidData& theData)
     // Prepare data for shrinking
     // ===========================
 
-    // Collect nodes to smooth, they are marked at the beginning of this method
+    // Collect nodes to smooth (they are marked at the beginning of this method)
     vector < const SMDS_MeshNode* > smoothNodes;
     {
       SMDS_NodeIteratorPtr nIt = smDS->GetNodes();
@@ -9816,8 +9863,13 @@ bool _ViscousBuilder::shrink(_SolidData& theData)
         if ( data._noShrinkShapes.count( subID ))
           continue;
         _EdgesOnShape* eos = data.GetShapeEdges( subID );
-        if ( !eos || eos->_sWOL.IsNull() ) continue;
-
+        if ( !eos || eos->_sWOL.IsNull() )
+          if ( data2 ) // check in adjacent SOLID
+          {
+            eos = data2->GetShapeEdges( subID );
+            if ( !eos || eos->_sWOL.IsNull() )
+              continue;
+          }
         subEOS.push_back( eos );
 
         for ( size_t i = 0; i < eos->_edges.size(); ++i )
@@ -9877,7 +9929,7 @@ bool _ViscousBuilder::shrink(_SolidData& theData)
       {
         const SMDS_MeshNode* n = smoothNodes[i];
         nodesToSmooth[ i ]._node = n;
-        // src nodes must be replaced by tgt nodes to have tgt nodes in _simplices
+        // src nodes must be already replaced by tgt nodes to have tgt nodes in _simplices
         _Simplex::GetSimplices( n, nodesToSmooth[ i ]._simplices, ignoreShapes, 0, sortSimplices);
         // fix up incorrect uv of nodes on the FACE
         helper.GetNodeUV( F, n, 0, &isOkUV);
@@ -10246,6 +10298,8 @@ bool _ViscousBuilder::shrink(_SolidData& theData)
 
     // Set an event listener to clear FACE sub-mesh together with SOLID sub-mesh
     VISCOUS_3D::ToClearSubWithMain( sm, data._solid );
+    if ( data2 )
+      VISCOUS_3D::ToClearSubWithMain( sm, data2->_solid );
 
   } // loop on FACES to srink mesh on
 
@@ -10324,7 +10378,7 @@ bool _ViscousBuilder::prepareEdgeToShrink( _LayerEdge&            edge,
 
     double uSrc = helper.GetNodeU( E, srcNode, n2 );
     double uTgt = helper.GetNodeU( E, tgtNode, srcNode );
-    double u2   = helper.GetNodeU( E, n2, srcNode );
+    double u2   = helper.GetNodeU( E, n2,      srcNode );
 
     //edge._pos.clear();
 
@@ -10859,7 +10913,7 @@ void _Shrinker1D::AddEdge( const _LayerEdge*   e,
     _done = false;
   }
   // check _LayerEdge
-  if ( e == _edges[0] || e == _edges[1] )
+  if ( e == _edges[0] || e == _edges[1] || e->_nodes.size() < 2 )
     return;
   if ( eos.SWOLType() != TopAbs_EDGE )
     throw SALOME_Exception(LOCALIZED("Wrong _LayerEdge is added"));
@@ -11060,7 +11114,8 @@ bool _ViscousBuilder::addBoundaryElements(_SolidData& data)
     for ( int iE = 1; iE <= geomEdges.Extent(); ++iE )
     {
       const TopoDS_Edge& E = TopoDS::Edge( geomEdges(iE));
-      if ( data._noShrinkShapes.count( getMeshDS()->ShapeToIndex( E )))
+      const TGeomID edgeID = getMeshDS()->ShapeToIndex( E );
+      if ( data._noShrinkShapes.count( edgeID ))
         continue;
 
       // Get _LayerEdge's based on E
@@ -11109,10 +11164,9 @@ bool _ViscousBuilder::addBoundaryElements(_SolidData& data)
       // Find out orientation and type of face to create
 
       bool reverse = false, isOnFace;
-      
-      map< TGeomID, TopoDS_Shape >::iterator e2f =
-        data._shrinkShape2Shape.find( getMeshDS()->ShapeToIndex( E ));
       TopoDS_Shape F;
+
+      map< TGeomID, TopoDS_Shape >::iterator e2f = data._shrinkShape2Shape.find( edgeID );
       if (( isOnFace = ( e2f != data._shrinkShape2Shape.end() )))
       {
         F = e2f->second.Oriented( TopAbs_FORWARD );
@@ -11138,18 +11192,44 @@ bool _ViscousBuilder::addBoundaryElements(_SolidData& data)
       if ( !sm )
         return error("error in addBoundaryElements()", data._index);
 
+      // Find a proxy sub-mesh of the FACE of an adjacent SOLID, which will use the new boundary
+      // faces for 3D meshing (PAL23414)
+      SMESHDS_SubMesh* adjSM = 0;
+      if ( isOnFace )
+      {
+        const TGeomID   faceID = sm->GetID();
+        PShapeIteratorPtr soIt = helper.GetAncestors( F, *_mesh, TopAbs_SOLID );
+        while ( const TopoDS_Shape* solid = soIt->next() )
+          if ( !solid->IsSame( data._solid ))
+          {
+            size_t iData = _solids.FindIndex( *solid ) - 1;
+            if ( iData < _sdVec.size() &&
+                 _sdVec[ iData ]._ignoreFaceIds.count( faceID ) &&
+                 _sdVec[ iData ]._shrinkShape2Shape.count( edgeID ) == 0 )
+            {
+              SMESH_ProxyMesh::SubMesh* proxySub =
+                _sdVec[ iData ]._proxyMesh->getFaceSubM( TopoDS::Face( F ), /*create=*/false);
+              if ( proxySub && proxySub->NbElements() > 0 )
+                adjSM = proxySub;
+            }
+          }
+      }
+
       // Make faces
       const int dj1 = reverse ? 0 : 1;
       const int dj2 = reverse ? 1 : 0;
+      vector< const SMDS_MeshElement*> ff; // new faces row
+      SMESHDS_Mesh* m = getMeshDS();
       for ( size_t j = 1; j < ledges.size(); ++j )
       {
         vector< const SMDS_MeshNode*>&  nn1 = ledges[j-dj1]->_nodes;
         vector< const SMDS_MeshNode*>&  nn2 = ledges[j-dj2]->_nodes;
+        ff.resize( std::max( nn1.size(), nn2.size() ), NULL );
         if ( nn1.size() == nn2.size() )
         {
           if ( isOnFace )
             for ( size_t z = 1; z < nn1.size(); ++z )
-              sm->AddElement( getMeshDS()->AddFace( nn1[z-1], nn2[z-1], nn2[z], nn1[z] ));
+              sm->AddElement( ff[z-1] = m->AddFace( nn1[z-1], nn2[z-1], nn2[z], nn1[z] ));
           else
             for ( size_t z = 1; z < nn1.size(); ++z )
               sm->AddElement( new SMDS_FaceOfNodes( nn1[z-1], nn2[z-1], nn2[z], nn1[z] ));
@@ -11158,7 +11238,7 @@ bool _ViscousBuilder::addBoundaryElements(_SolidData& data)
         {
           if ( isOnFace )
             for ( size_t z = 1; z < nn2.size(); ++z )
-              sm->AddElement( getMeshDS()->AddFace( nn1[0], nn2[z-1], nn2[z] ));
+              sm->AddElement( ff[z-1] = m->AddFace( nn1[0], nn2[z-1], nn2[z] ));
           else
             for ( size_t z = 1; z < nn2.size(); ++z )
               sm->AddElement( new SMDS_FaceOfNodes( nn1[0], nn2[z-1], nn2[z] ));
@@ -11167,11 +11247,19 @@ bool _ViscousBuilder::addBoundaryElements(_SolidData& data)
         {
           if ( isOnFace )
             for ( size_t z = 1; z < nn1.size(); ++z )
-              sm->AddElement( getMeshDS()->AddFace( nn1[z-1], nn2[0], nn1[z] ));
+              sm->AddElement( ff[z-1] = m->AddFace( nn1[z-1], nn2[0], nn1[z] ));
           else
             for ( size_t z = 1; z < nn1.size(); ++z )
               sm->AddElement( new SMDS_FaceOfNodes( nn1[z-1], nn2[0], nn2[z] ));
         }
+
+        if ( adjSM ) // add faces to a proxy SM of the adjacent SOLID
+        {
+          for ( size_t z = 0; z < ff.size(); ++z )
+            if ( ff[ z ])
+              adjSM->AddElement( ff[ z ]);
+          ff.clear();
+        }
       }
 
       // Make edges
@@ -11182,7 +11270,7 @@ bool _ViscousBuilder::addBoundaryElements(_SolidData& data)
         if ( eos && eos->SWOLType() == TopAbs_EDGE )
         {
           vector< const SMDS_MeshNode*>&  nn = edge->_nodes;
-          if ( nn.size() < 2 || nn[1]->GetInverseElementIterator( SMDSAbs_Edge )->more() )
+          if ( nn.size() < 2 || nn[1]->NbInverseElements( SMDSAbs_Edge ) >= 2 )
             continue;
           helper.SetSubShape( eos->_sWOL );
           helper.SetElementsOnShape( true );