Salome HOME
Fix SALOME_TESTS/Grids/smesh/bugs_07/H5 as the geometry changed
[modules/smesh.git] / src / SMESH / SMESH_MesherHelper.cxx
index 913362e11fb3ac64601543b9a2e45161bf9f4e77..a1f9b3963a8c7be2ceabac23936ac6a3937e3ce5 100644 (file)
 #include "SMDS_FacePosition.hxx" 
 #include "SMDS_IteratorOnIterators.hxx"
 #include "SMDS_VolumeTool.hxx"
+#include "SMESHDS_Mesh.hxx"
 #include "SMESH_Block.hxx"
 #include "SMESH_HypoFilter.hxx"
+#include "SMESH_Mesh.hxx"
 #include "SMESH_MeshAlgos.hxx"
 #include "SMESH_ProxyMesh.hxx"
 #include "SMESH_subMesh.hxx"
@@ -113,6 +115,28 @@ SMESH_MesherHelper::~SMESH_MesherHelper()
   }
 }
 
+//================================================================================
+/*!
+ * \brief Return SMESH_Gen
+ */
+//================================================================================
+
+SMESH_Gen* SMESH_MesherHelper::GetGen() const
+{
+  return GetMesh()->GetGen();
+}
+
+//================================================================================
+/*!
+ * \brief Return mesh DS
+ */
+//================================================================================
+
+SMESHDS_Mesh* SMESH_MesherHelper::GetMeshDS() const
+{
+  return GetMesh()->GetMeshDS();
+}
+
 //=======================================================================
 //function : IsQuadraticSubMesh
 //purpose  : Check submesh for given shape: if all elements on this shape 
@@ -280,6 +304,7 @@ void SMESH_MesherHelper::SetSubShape(const TopoDS_Shape& aSh)
             double u2 = uv1.Coord(1);
             myPar1[0] = Min( u1, u2 );
             myPar2[0] = Max( u1, u2 );
+            myParIndex |= U_periodic;
           }
           else
           {
@@ -289,6 +314,7 @@ void SMESH_MesherHelper::SetSubShape(const TopoDS_Shape& aSh)
             double v2 = uv1.Coord(2);
             myPar1[1] = Min( v1, v2 );
             myPar2[1] = Max( v1, v2 );
+            myParIndex |= V_periodic;
           }
         }
         else //if ( !isSeam )
@@ -347,6 +373,16 @@ void SMESH_MesherHelper::SetSubShape(const TopoDS_Shape& aSh)
   }
 }
 
+//=======================================================================
+//function : ShapeToIndex
+//purpose  : Convert a shape to its index in the SMESHDS_Mesh
+//=======================================================================
+
+int SMESH_MesherHelper::ShapeToIndex( const TopoDS_Shape& S ) const
+{
+  return GetMeshDS()->ShapeToIndex( S );
+}
+
 //=======================================================================
 //function : GetNodeUVneedInFaceNode
 //purpose  : Check if inFaceNode argument is necessary for call GetNodeUV(F,..)
@@ -1613,7 +1649,7 @@ const SMDS_MeshNode* SMESH_MesherHelper::GetMediumNode(const SMDS_MeshNode* n1,
 
   // get type of shape for the new medium node
   int faceID = -1, edgeID = -1;
-  TopoDS_Edge E; double u [2];
+  TopoDS_Edge E; double u [2] = {0.,0.};
   TopoDS_Face F; gp_XY  uv[2];
   bool uvOK[2] = { true, true };
   const bool useCurSubShape = ( !myShape.IsNull() && myShape.ShapeType() == TopAbs_EDGE );
@@ -2532,6 +2568,7 @@ bool SMESH_MesherHelper::LoadNodeColumns(TParam2ColumnMap &            theParam2
     }
 
     // get nodes on theBaseEdge sorted by param on edge and initialize theParam2ColumnMap with them
+    const SMDS_MeshNode* prevEndNodes[2] = { 0, 0 };
     edge = theBaseSide.begin();
     for ( int iE = 0; edge != theBaseSide.end(); ++edge, ++iE )
     {
@@ -2599,11 +2636,16 @@ bool SMESH_MesherHelper::LoadNodeColumns(TParam2ColumnMap &            theParam2
       const double prevPar = theParam2ColumnMap.empty() ? 0 : theParam2ColumnMap.rbegin()->first;
       for ( u_n = sortedBaseNN.begin(); u_n != sortedBaseNN.end(); u_n++ )
       {
+        if ( u_n->second == prevEndNodes[0] ||
+             u_n->second == prevEndNodes[1] )
+          continue;
         double par = prevPar + coeff * ( u_n->first - f );
         TParam2ColumnMap::iterator u2nn =
           theParam2ColumnMap.insert( theParam2ColumnMap.end(), make_pair( par, TNodeColumn()));
         u2nn->second.push_back( u_n->second );
       }
+      prevEndNodes[0] = sortedBaseNN.begin()->second;
+      prevEndNodes[1] = sortedBaseNN.rbegin()->second;
     }
     if ( theParam2ColumnMap.size() < 2 )
       return false;
@@ -2894,7 +2936,7 @@ bool SMESH_MesherHelper::IsReversedSubMesh (const TopoDS_Face& theFace)
         TopoDS_Shape s0 = GetSubShapeByNode( nn[0], GetMeshDS() );
         TopoDS_Shape s1 = GetSubShapeByNode( nn[1], GetMeshDS() );
         TopoDS_Shape  E = GetCommonAncestor( s0, s1, *myMesh, TopAbs_EDGE );
-        if ( !E.IsNull() && !s0.IsSame( s1 ))
+        if ( !E.IsNull() && !s0.IsSame( s1 ) && E.Orientation() != TopAbs_INTERNAL )
         {
           // is E seam edge?
           int nb = 0;
@@ -2908,10 +2950,16 @@ bool SMESH_MesherHelper::IsReversedSubMesh (const TopoDS_Face& theFace)
             bool ok = true;
             double u0 = GetNodeU( TopoDS::Edge( E ), nn[0], nn[1], &ok );
             double u1 = GetNodeU( TopoDS::Edge( E ), nn[1], nn[0], &ok );
-            // check that the 2 nodes are connected with a segment (IPAL53055)
-            if ( SMESHDS_SubMesh* sm = GetMeshDS()->MeshElements( E ))
-              if ( sm->NbElements() > 0 && !GetMeshDS()->FindEdge( nn[0], nn[1] ))
-                ok = false;
+            if ( ok )
+            {
+              // check that the 2 nodes are connected with a segment (IPAL53055)
+              ok = false;
+              const SMDS_MeshElement* seg;
+              if ( SMESHDS_SubMesh* sm = GetMeshDS()->MeshElements( E ))
+                if (( sm->NbElements() > 0 ) &&
+                    ( seg = GetMeshDS()->FindEdge( nn[0], nn[1] )))
+                  ok = sm->Contains( seg );
+            }
             if ( ok )
             {
               isReversed = ( u0 > u1 );
@@ -3360,11 +3408,16 @@ namespace {
     TopTools_ListIteratorOfListOfShape _ancIter;
     TopAbs_ShapeEnum                   _type;
     TopTools_MapOfShape                _encountered;
-    TAncestorsIterator( const TopTools_ListOfShape& ancestors, TopAbs_ShapeEnum type)
+    TopTools_IndexedMapOfShape         _allowed;
+    TAncestorsIterator( const TopTools_ListOfShape& ancestors,
+                        TopAbs_ShapeEnum            type,
+                        const TopoDS_Shape*         container/* = 0*/)
       : _ancIter( ancestors ), _type( type )
     {
+      if ( container && !container->IsNull() )
+        TopExp::MapShapes( *container, type, _allowed);
       if ( _ancIter.More() ) {
-        if ( _ancIter.Value().ShapeType() != _type ) next();
+        if ( !isCurrentAllowed() ) next();
         else _encountered.Add( _ancIter.Value() );
       }
     }
@@ -3377,25 +3430,32 @@ namespace {
       const TopoDS_Shape* s = _ancIter.More() ? & _ancIter.Value() : 0;
       if ( _ancIter.More() )
         for ( _ancIter.Next();  _ancIter.More(); _ancIter.Next())
-          if ( _ancIter.Value().ShapeType() == _type && _encountered.Add( _ancIter.Value() ))
+          if ( isCurrentAllowed() && _encountered.Add( _ancIter.Value() ))
             break;
       return s;
     }
+    bool isCurrentAllowed()
+    {
+      return (( _ancIter.Value().ShapeType() == _type ) &&
+              ( _allowed.IsEmpty() || _allowed.Contains( _ancIter.Value() )));
+    }
   };
 
 } // namespace
 
 //=======================================================================
 /*!
- * \brief Return iterator on ancestors of the given type
+ * \brief Return iterator on ancestors of the given type, included into a container shape
  */
 //=======================================================================
 
 PShapeIteratorPtr SMESH_MesherHelper::GetAncestors(const TopoDS_Shape& shape,
                                                    const SMESH_Mesh&   mesh,
-                                                   TopAbs_ShapeEnum    ancestorType)
+                                                   TopAbs_ShapeEnum    ancestorType,
+                                                   const TopoDS_Shape* container)
 {
-  return PShapeIteratorPtr( new TAncestorsIterator( mesh.GetAncestors(shape), ancestorType));
+  return PShapeIteratorPtr
+    ( new TAncestorsIterator( mesh.GetAncestors(shape), ancestorType, container));
 }
 
 //=======================================================================
@@ -4593,7 +4653,7 @@ namespace { // Structures used by FixQuadraticElements()
       SMDS_ElemIteratorPtr faceIter( new TIterOnIter( faceIterVec ));
 
       // a seacher to check if a volume is close to a concave face
-      std::auto_ptr< SMESH_ElementSearcher > faceSearcher
+      SMESHUtils::Deleter< SMESH_ElementSearcher > faceSearcher
         ( SMESH_MeshAlgos::GetElementSearcher( *theHelper.GetMeshDS(), faceIter ));
 
       // classifier
@@ -4695,7 +4755,7 @@ namespace { // Structures used by FixQuadraticElements()
                     gp_Pnt pMedium = SMESH_TNodeXYZ( linkIt->second );
                     double hMedium = faceNorm * gp_Vec( pOnFace0, pMedium ).XYZ();
                     double hVol    = faceNorm * gp_Vec( pOnFace0, pInSolid ).XYZ();
-                    isDistorted = ( Abs( hMedium ) > Abs( hVol * 0.5 ));
+                    isDistorted = ( Abs( hMedium ) > Abs( hVol * 0.75 ));
                   }
                 }
               }
@@ -4776,7 +4836,8 @@ void SMESH_MesherHelper::FixQuadraticElements(SMESH_ComputeErrorPtr& compError,
     }
     // fix nodes on geom faces
 #ifdef _DEBUG_
-    int nbfaces = faces.Extent(); /*avoid "unused varianbles": */ nbfaces++, nbfaces--; 
+    int nbfaces = nbSolids;
+    nbfaces = faces.Extent(); /*avoid "unused varianbles": */ nbfaces++, nbfaces--; 
 #endif
     for ( TopTools_MapIteratorOfMapOfShape fIt( faces ); fIt.More(); fIt.Next() ) {
       MSG("FIX FACE " << nbfaces-- << " #" << GetMeshDS()->ShapeToIndex(fIt.Key()));
@@ -5117,6 +5178,7 @@ void SMESH_MesherHelper::FixQuadraticElements(SMESH_ComputeErrorPtr& compError,
                      "uv2: "<<uv2.X()<<", "<<uv2.Y()<<" \t" <<
                      "uvOld: "<<oldUV.X()<<", "<<oldUV.Y()<<" \t" <<
                      "newUV: "<<newUV.X()<<", "<<newUV.Y()<<" \t");
+                uv0.SetX( uv2.X() ); // avoid warning: variable set but not used
               }
 #endif
               (*link1)->Move( move, /*sum=*/false, /*is2dFixed=*/true );
@@ -5135,7 +5197,6 @@ void SMESH_MesherHelper::FixQuadraticElements(SMESH_ComputeErrorPtr& compError,
   // -------------
 
   TIDSortedElemSet biQuadQuas, biQuadTris, triQuadHexa;
-  const SMDS_MeshElement *biQuadQua, *triQuadHex;
   const bool toFixCentralNodes = ( myMesh->NbBiQuadQuadrangles() +
                                    myMesh->NbBiQuadTriangles() +
                                    myMesh->NbTriQuadraticHexas() );
@@ -5166,7 +5227,6 @@ void SMESH_MesherHelper::FixQuadraticElements(SMESH_ComputeErrorPtr& compError,
       // collect bi-quadratic elements
       if ( toFixCentralNodes )
       {
-        biQuadQua = triQuadHex = 0;
         SMDS_ElemIteratorPtr eIt = pLink->_mediumNode->GetInverseElementIterator();
         while ( eIt->more() )
         {
@@ -5325,6 +5385,10 @@ void SMESH_MesherHelper::FixQuadraticElements(SMESH_ComputeErrorPtr& compError,
                              nCenterCoords.X(), nCenterCoords.Y(), nCenterCoords.Z());
     }
   }
+#ifdef _DEBUG_
+  // avoid warning: defined but not used operator<<()
+  SMESH_Comment() << *links.begin() << *faces.begin();
+#endif
 }
 
 //================================================================================