Salome HOME
Merge from V6_main 12/11/2012
[modules/smesh.git] / src / StdMeshers / StdMeshers_Import_1D2D.cxx
index 29ffeaf9acc6ea80509186140378a8c7f3374106..8180f0a027ff93a692f6e7155cd949443c4f0902 100644 (file)
@@ -1,23 +1,23 @@
-//  Copyright (C) 2007-2010  CEA/DEN, EDF R&D, OPEN CASCADE
+// Copyright (C) 2007-2012  CEA/DEN, EDF R&D, OPEN CASCADE
 //
-//  Copyright (C) 2003-2007  OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
-//  CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS
+// Copyright (C) 2003-2007  OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
+// CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS
 //
-//  This library is free software; you can redistribute it and/or
-//  modify it under the terms of the GNU Lesser General Public
-//  License as published by the Free Software Foundation; either
-//  version 2.1 of the License.
+// This library is free software; you can redistribute it and/or
+// modify it under the terms of the GNU Lesser General Public
+// License as published by the Free Software Foundation; either
+// version 2.1 of the License.
 //
-//  This library is distributed in the hope that it will be useful,
-//  but WITHOUT ANY WARRANTY; without even the implied warranty of
-//  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
-//  Lesser General Public License for more details.
+// This library is distributed in the hope that it will be useful,
+// but WITHOUT ANY WARRANTY; without even the implied warranty of
+// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+// Lesser General Public License for more details.
 //
-//  You should have received a copy of the GNU Lesser General Public
-//  License along with this library; if not, write to the Free Software
-//  Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
+// You should have received a copy of the GNU Lesser General Public
+// License along with this library; if not, write to the Free Software
+// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
 //
-//  See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
+// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
 //
 
 //  SMESH SMESH : implementaion of SMESH idl descriptions
 #include <TopoDS_Compound.hxx>
 #include <TopoDS_Edge.hxx>
 #include <TopoDS_Vertex.hxx>
+#include <Precision.hxx>
 
 #include <numeric>
 
 using namespace std;
 
+namespace
+{
+  double getMinElemSize2( const SMESHDS_GroupBase* srcGroup )
+  {
+    double minSize2 = 1e100;
+    SMDS_ElemIteratorPtr srcElems = srcGroup->GetElements();
+    while ( srcElems->more() ) // loop on group contents
+    {
+      const SMDS_MeshElement* face = srcElems->next();
+      int nbN = face->NbCornerNodes();
+
+      SMESH_TNodeXYZ prevN( face->GetNode( nbN-1 ));
+      for ( int i = 0; i < nbN; ++i )
+      {
+        SMESH_TNodeXYZ n( face->GetNode( i ) );
+        double size2 = ( n - prevN ).SquareModulus();
+        minSize2 = std::min( minSize2, size2 );
+        prevN = n;
+      }
+    }
+    return minSize2;
+  }
+}
+
 //=============================================================================
 /*!
  * Creates StdMeshers_Import_1D2D
@@ -70,7 +95,7 @@ StdMeshers_Import_1D2D::StdMeshers_Import_1D2D(int hypId, int studyId, SMESH_Gen
   _shapeType = (1 << TopAbs_FACE);
 
   _compatibleHypothesis.push_back("ImportSource2D");
-  _requireDescretBoundary = false;
+  _requireDiscreteBoundary = false;
 }
 
 //=============================================================================
@@ -143,6 +168,12 @@ bool StdMeshers_Import_1D2D::Compute(SMESH_Mesh & theMesh, const TopoDS_Shape &
   if ( srcGroups.empty() )
     return error("Invalid source groups");
 
+  bool allGroupsEmpty = true;
+  for ( size_t iG = 0; iG < srcGroups.size() && allGroupsEmpty; ++iG )
+    allGroupsEmpty = srcGroups[iG]->GetGroupDS()->IsEmpty();
+  if ( allGroupsEmpty )
+    return error("No faces in source groups");
+
   SMESH_MesherHelper helper(theMesh);
   helper.SetSubShape(theShape);
   SMESHDS_Mesh* tgtMesh = theMesh.GetMeshDS();
@@ -182,6 +213,7 @@ bool StdMeshers_Import_1D2D::Compute(SMESH_Mesh & theMesh, const TopoDS_Shape &
   // to count now many times a link between nodes encounters
   map<TLink, int> linkCount;
   map<TLink, int>::iterator link2Nb;
+  double minGroupTol = Precision::Infinite();
 
   // =========================
   // Import faces from groups
@@ -190,7 +222,7 @@ bool StdMeshers_Import_1D2D::Compute(SMESH_Mesh & theMesh, const TopoDS_Shape &
   StdMeshers_Import_1D::TNodeNodeMap* n2n;
   StdMeshers_Import_1D::TElemElemMap* e2e;
   vector<const SMDS_MeshNode*> newNodes;
-  for ( int iG = 0; iG < srcGroups.size(); ++iG )
+  for ( size_t iG = 0; iG < srcGroups.size(); ++iG )
   {
     const SMESHDS_GroupBase* srcGroup = srcGroups[iG]->GetGroupDS();
 
@@ -199,9 +231,12 @@ bool StdMeshers_Import_1D2D::Compute(SMESH_Mesh & theMesh, const TopoDS_Shape &
     if ( !srcMesh ) continue;
     StdMeshers_Import_1D::getMaps( srcMesh, &theMesh, n2n, e2e );
 
+    const double groupTol = 0.5 * sqrt( getMinElemSize2( srcGroup ));
+    minGroupTol = std::min( groupTol, minGroupTol );
+
     SMDS_ElemIteratorPtr srcElems = srcGroup->GetElements();
     SMDS_MeshNode *tmpNode = helper.AddNode(0,0,0);
-    gp_XY uv;
+    gp_XY uv( Precision::Infinite(), Precision::Infinite() );
     while ( srcElems->more() ) // loop on group contents
     {
       const SMDS_MeshElement* face = srcElems->next();
@@ -210,9 +245,9 @@ bool StdMeshers_Import_1D2D::Compute(SMESH_Mesh & theMesh, const TopoDS_Shape &
       newNodes.back() = 0;
       int nbCreatedNodes = 0;
       SMDS_MeshElement::iterator node = face->begin_nodes();
-      for ( unsigned i = 0; i < newNodes.size(); ++i, ++node )
+      for ( size_t i = 0; i < newNodes.size(); ++i, ++node )
       {
-        TNodeNodeMap::iterator n2nIt = n2n->insert( make_pair( *node, (SMDS_MeshNode*)0 )).first;
+        StdMeshers_Import_1D::TNodeNodeMap::iterator n2nIt = n2n->insert( make_pair( *node, (SMDS_MeshNode*)0 )).first;
         if ( n2nIt->second )
         {
           if ( !subShapeIDs.count( n2nIt->second->getshapeId() ))
@@ -222,7 +257,7 @@ bool StdMeshers_Import_1D2D::Compute(SMESH_Mesh & theMesh, const TopoDS_Shape &
         {
           // find an existing vertex node
           for ( vNIt = vertexNodes.begin(); vNIt != vertexNodes.end(); ++vNIt)
-            if ( vNIt->SquareDistance( *node ) < 10 * faceTol * faceTol)
+            if ( vNIt->SquareDistance( *node ) < groupTol * groupTol)
             {
               (*n2nIt).second = vNIt->_node;
               vertexNodes.erase( vNIt );
@@ -233,7 +268,8 @@ bool StdMeshers_Import_1D2D::Compute(SMESH_Mesh & theMesh, const TopoDS_Shape &
         {
           // find out if node lies on theShape
           tmpNode->setXYZ( (*node)->X(), (*node)->Y(), (*node)->Z());
-          if ( helper.CheckNodeUV( geomFace, tmpNode, uv, 10 * faceTol, /*force=*/true ))
+          uv.SetCoord( Precision::Infinite(), Precision::Infinite() );
+          if ( helper.CheckNodeUV( geomFace, tmpNode, uv, groupTol, /*force=*/true ))
           {
             SMDS_MeshNode* newNode = tgtMesh->AddNode( (*node)->X(), (*node)->Y(), (*node)->Z());
             n2nIt->second = newNode;
@@ -313,6 +349,13 @@ bool StdMeshers_Import_1D2D::Compute(SMESH_Mesh & theMesh, const TopoDS_Shape &
           medium = newNodes[i+nbNodes];
         link2Nb = linkCount.insert( make_pair( TLink( n1, n2, medium ), 0)).first;
         ++link2Nb->second;
+        // if ( link2Nb->second == 1 )
+        // {
+        //   // measure link length
+        //   double len2 = SMESH_TNodeXYZ( n1 ).SquareDistance( n2 );
+        //   if ( len2 < minGroupTol )
+        //     minGroupTol = len2;
+        // }
       }
     }
     helper.GetMeshDS()->RemoveNode(tmpNode);
@@ -328,8 +371,13 @@ bool StdMeshers_Import_1D2D::Compute(SMESH_Mesh & theMesh, const TopoDS_Shape &
     if ( subShapeIDs.insert( tgtMesh->ShapeToIndex( exp.Current() )).second )
       edges.push_back( TopoDS::Edge( exp.Current() ));
 
+  // use large tolerance for projection of nodes to edges because of
+  // BLSURF mesher specifics (issue 0020918, Study2.hdf)
+  const double projTol = minGroupTol;
+
   bool isFaceMeshed = false;
-  if ( SMESHDS_SubMesh* tgtSM = tgtMesh->MeshElements( theShape ))
+  SMESHDS_SubMesh* tgtFaceSM = tgtMesh->MeshElements( theShape );
+  if ( tgtFaceSM )
   {
     // the imported mesh is valid if all external links (encountered once)
     // lie on geom edges
@@ -341,7 +389,7 @@ bool StdMeshers_Import_1D2D::Compute(SMESH_Mesh & theMesh, const TopoDS_Shape &
       int nbFaces = link2Nb->second;
       if ( nbFaces == 1 )
       {
-        // check if a not shared link lie on face boundary
+        // check if a not shared link lies on face boundary
         bool nodesOnBoundary = true;
         list< TopoDS_Shape > bndShapes;
         for ( int is1stN = 0; is1stN < 2 && nodesOnBoundary; ++is1stN )
@@ -349,14 +397,14 @@ bool StdMeshers_Import_1D2D::Compute(SMESH_Mesh & theMesh, const TopoDS_Shape &
           const SMDS_MeshNode* n = is1stN ? link.node1() : link.node2();
           if ( !subShapeIDs.count( n->getshapeId() ))
           {
-            for ( unsigned iE = 0; iE < edges.size(); ++iE )
-              if ( helper.CheckNodeU( edges[iE], n, u, 10 * faceTol, /*force=*/true ))
+            for ( size_t iE = 0; iE < edges.size(); ++iE )
+              if ( helper.CheckNodeU( edges[iE], n, u=0, projTol, /*force=*/true ))
               {
                 BRep_Tool::Range(edges[iE],f,l);
                 if ( Abs(u-f) < 2 * faceTol || Abs(u-l) < 2 * faceTol )
                   // duplicated node on vertex
                   return error("Source elements overlap one another");
-                tgtSM->RemoveNode( n, /*isNodeDeleted=*/false );
+                tgtFaceSM->RemoveNode( n, /*isNodeDeleted=*/false );
                 tgtMesh->SetNodeOnEdge( (SMDS_MeshNode*)n, edges[iE], u );
                 break;
               }
@@ -372,25 +420,27 @@ bool StdMeshers_Import_1D2D::Compute(SMESH_Mesh & theMesh, const TopoDS_Shape &
           }
         }
         if ( !nodesOnBoundary )
-          break; // error: free internal link
+        {
+          error("free internal link"); // just for an easier debug
+          break;
+        }
         if ( bndShapes.front().ShapeType() == TopAbs_EDGE &&
              bndShapes.front() != bndShapes.back() )
-          break; // error: link nodes on different geom edges
+          // link nodes on different geom edges
+          return error(COMPERR_BAD_INPUT_MESH, "Source nodes mismatch target vertices");
 
         // find geom edge the link is on
         if ( bndShapes.back().ShapeType() != TopAbs_EDGE )
         {
           // find geom edge by two vertices
-          TopoDS_Shape geomEdge;
-          PShapeIteratorPtr edgeIt = helper.GetAncestors( bndShapes.back(), theMesh, TopAbs_EDGE );
-          while ( edgeIt->more() )
+          TopoDS_Shape geomEdge = helper.GetCommonAncestor( bndShapes.back(),
+                                                            bndShapes.front(),
+                                                            theMesh, TopAbs_EDGE );
+          if ( geomEdge.IsNull() )
           {
-            geomEdge = *(edgeIt->next());
-            if ( !helper.IsSubShape( bndShapes.front(), geomEdge ))
-              geomEdge.Nullify();
+            error("free internal link");
+            break; // vertices belong to different edges
           }
-          if ( geomEdge.IsNull() )
-            break; // vertices belong to different edges -> error: free internal link
           bndShapes.push_back( geomEdge );
         }
 
@@ -407,8 +457,8 @@ bool StdMeshers_Import_1D2D::Compute(SMESH_Mesh & theMesh, const TopoDS_Shape &
           edge = tgtMesh->AddEdge( newNodes[0], newNodes[1], newNodes[2] );
 
           TopoDS_Edge geomEdge = TopoDS::Edge(bndShapes.back());
-          helper.CheckNodeU( geomEdge, link._medium, u, 10*faceTol, /*force=*/true );
-          tgtSM->RemoveNode( link._medium, /*isNodeDeleted=*/false );
+          helper.CheckNodeU( geomEdge, link._medium, u, projTol, /*force=*/true );
+          tgtFaceSM->RemoveNode( link._medium, /*isNodeDeleted=*/false );
           tgtMesh->SetNodeOnEdge( (SMDS_MeshNode*)link._medium, geomEdge, u );
         }
         else
@@ -422,7 +472,7 @@ bool StdMeshers_Import_1D2D::Compute(SMESH_Mesh & theMesh, const TopoDS_Shape &
       }
       else if ( nbFaces > 2 )
       {
-        return error( "Non-manifold source mesh");
+        return error( COMPERR_BAD_INPUT_MESH, "Non-manifold source mesh");
       }
     }
     isFaceMeshed = ( link2Nb == linkCount.end() && !linkCount.empty());
@@ -448,23 +498,105 @@ bool StdMeshers_Import_1D2D::Compute(SMESH_Mesh & theMesh, const TopoDS_Shape &
         if ( nbEdges < 2 )
           return false; // weird
         if ( nbEdges > 2 )
-          return error( "Source elements overlap one another");
+          return error( COMPERR_BAD_INPUT_MESH, "Source elements overlap one another");
       }
     }
   }
   if ( !isFaceMeshed )
-    return error( "Source elements don't cover totally the geometrical face" );
+    return error( COMPERR_BAD_INPUT_MESH,
+                  "Source elements don't cover totally the geometrical face" );
+
+  if ( helper.HasSeam() )
+  {
+    // links on seam edges are shared by two faces, so no edges were created on them
+    // by the previous detection of 2D mesh boundary
+    for ( size_t iE = 0; iE < edges.size(); ++iE )
+    {
+      if ( !helper.IsRealSeam( edges[iE] )) continue;
+      const TopoDS_Edge& seamEdge = edges[iE];
+      // to find nodes lying on the seamEdge we check nodes of mesh faces sharing a node on one
+      // of its vertices; after finding another node on seamEdge we continue the same way
+      // until finding all nodes.
+      TopoDS_Vertex      seamVertex = helper.IthVertex( 0, seamEdge );
+      const SMDS_MeshNode* vertNode = SMESH_Algo::VertexNode( seamVertex, tgtMesh );
+      set< const SMDS_MeshNode* > checkedNodes; checkedNodes.insert( vertNode );
+      set< const SMDS_MeshElement* > checkedFaces;
+      // as a face can have more than one node on the seamEdge, there is a difficulty in selecting
+      // one of those nodes to treat next; so we simply find all nodes on the seamEdge and
+      // then sort them by U on edge
+      typedef list< pair< double, const SMDS_MeshNode* > > TUNodeList;
+      TUNodeList nodesOnSeam;
+      double u = helper.GetNodeU( seamEdge, vertNode );
+      nodesOnSeam.push_back( make_pair( u, vertNode ));
+      TUNodeList::iterator u2nIt = nodesOnSeam.begin();
+      for ( ; u2nIt != nodesOnSeam.end(); ++u2nIt )
+      {
+        const SMDS_MeshNode* startNode = (*u2nIt).second;
+        SMDS_ElemIteratorPtr faceIt = startNode->GetInverseElementIterator( SMDSAbs_Face );
+        while ( faceIt->more() )
+        {
+          const SMDS_MeshElement* face = faceIt->next();
+          if ( !checkedFaces.insert( face ).second ) continue;
+          for ( int i = 0, nbNodes = face->NbCornerNodes(); i < nbNodes; ++i )
+          {
+            const SMDS_MeshNode* n = face->GetNode( i );
+            if ( n == startNode || !checkedNodes.insert( n ).second ) continue;
+            if ( helper.CheckNodeU( seamEdge, n, u=0, projTol, /*force=*/true ))
+              nodesOnSeam.push_back( make_pair( u, n ));
+          }
+        }
+      }
+      // sort the found nodes by U on the seamEdge; most probably they are in a good order,
+      // so we can use the hint to spead-up map filling
+      map< double, const SMDS_MeshNode* > u2nodeMap;
+      for ( u2nIt = nodesOnSeam.begin(); u2nIt != nodesOnSeam.end(); ++u2nIt )
+        u2nodeMap.insert( u2nodeMap.end(), *u2nIt );
+
+      // create edges
+      {
+        SMESH_MesherHelper seamHelper( theMesh );
+        seamHelper.SetSubShape( edges[ iE ]);
+        seamHelper.SetElementsOnShape( true );
+
+        if ( (*checkedFaces.begin())->IsQuadratic() )
+          for ( set< const SMDS_MeshElement* >::iterator fIt = checkedFaces.begin();
+                fIt != checkedFaces.end(); ++fIt )
+            seamHelper.AddTLinks( static_cast<const SMDS_MeshFace*>( *fIt ));
+
+        map< double, const SMDS_MeshNode* >::iterator n1, n2, u2nEnd = u2nodeMap.end();
+        for ( n2 = u2nodeMap.begin(), n1 = n2++; n2 != u2nEnd; ++n1, ++n2 )
+        {
+          const SMDS_MeshNode* node1 = n1->second;
+          const SMDS_MeshNode* node2 = n2->second;
+          seamHelper.AddEdge( node1, node2 );
+          if ( node2->getshapeId() == helper.GetSubShapeID() )
+          {
+            tgtFaceSM->RemoveNode( node2, /*isNodeDeleted=*/false );
+            tgtMesh->SetNodeOnEdge( const_cast<SMDS_MeshNode*>( node2 ), seamEdge, n2->first );
+          }
+        }
+      }
+    } // loop on edges to find seam ones
+  } // if ( helper.HasSeam() )
 
   // notify sub-meshes of edges on computation
-  for ( unsigned iE = 0; iE < edges.size(); ++iE )
-    theMesh.GetSubMesh( edges[iE] )->ComputeStateEngine(SMESH_subMesh::CHECK_COMPUTE_STATE);
+  for ( size_t iE = 0; iE < edges.size(); ++iE )
+  {
+    SMESH_subMesh * sm = theMesh.GetSubMesh( edges[iE] );
+    if ( BRep_Tool::Degenerated( edges[iE] ))
+      sm->SetIsAlwaysComputed( true );
+    sm->ComputeStateEngine(SMESH_subMesh::CHECK_COMPUTE_STATE);
+    if ( sm->GetComputeState() != SMESH_subMesh::COMPUTE_OK )
+      return error(SMESH_Comment("Failed to create segments on the edge ")
+                   << tgtMesh->ShapeToIndex( edges[iE ]));
+  }
 
   // ============
   // Copy meshes
   // ============
 
   vector<SMESH_Mesh*> srcMeshes = _sourceHyp->GetSourceMeshes();
-  for ( unsigned i = 0; i < srcMeshes.size(); ++i )
+  for ( size_t i = 0; i < srcMeshes.size(); ++i )
     StdMeshers_Import_1D::importMesh( srcMeshes[i], theMesh, _sourceHyp, theShape );
 
   return true;
@@ -535,7 +667,6 @@ bool StdMeshers_Import_1D2D::Evaluate(SMESH_Mesh &         theMesh,
     helper.SetSubShape(theShape);
 
     const TopoDS_Face& geomFace = TopoDS::Face( theShape );
-    const double faceTol = helper.MaxTolerance( geomFace );
 
     // take into account nodes on vertices
     TopExp_Explorer exp( theShape, TopAbs_VERTEX );
@@ -550,9 +681,12 @@ bool StdMeshers_Import_1D2D::Evaluate(SMESH_Mesh &         theMesh,
     // count faces and nodes imported from groups
     set<const SMDS_MeshNode* > allNodes;
     gp_XY uv;
+    double minGroupTol = 1e100;
     for ( int iG = 0; iG < srcGroups.size(); ++iG )
     {
       const SMESHDS_GroupBase* srcGroup = srcGroups[iG]->GetGroupDS();
+      const double groupTol = 0.5 * sqrt( getMinElemSize2( srcGroup ));
+      minGroupTol = std::min( groupTol, minGroupTol );
       SMDS_ElemIteratorPtr srcElems = srcGroup->GetElements();
       SMDS_MeshNode *tmpNode =helper.AddNode(0,0,0);
       while ( srcElems->more() ) // loop on group contents
@@ -563,7 +697,7 @@ bool StdMeshers_Import_1D2D::Evaluate(SMESH_Mesh &         theMesh,
         gp_XYZ gc(0,0,0);
         gc = accumulate( TXyzIterator(face->nodesIterator()), TXyzIterator(), gc)/face->NbNodes();
         tmpNode->setXYZ( gc.X(), gc.Y(), gc.Z());
-        if ( helper.CheckNodeUV( geomFace, tmpNode, uv, 10 * faceTol, /*force=*/true ))
+        if ( helper.CheckNodeUV( geomFace, tmpNode, uv, groupTol, /*force=*/true ))
         {
           ++aVec[ face->GetEntityType() ];
 
@@ -609,8 +743,8 @@ bool StdMeshers_Import_1D2D::Evaluate(SMESH_Mesh &         theMesh,
           bool eraseLink = ( nbFacesOfLink != 1 );
           if ( nbFacesOfLink == 1 )
           {
-            if ( helper.CheckNodeU( geomEdge, link.node1(), u, 10*faceTol, /*force=*/true )&&
-                 helper.CheckNodeU( geomEdge, link.node2(), u, 10*faceTol, /*force=*/true ))
+            if ( helper.CheckNodeU( geomEdge, link.node1(), u, minGroupTol, /*force=*/true )&&
+                 helper.CheckNodeU( geomEdge, link.node2(), u, minGroupTol, /*force=*/true ))
             {
               bool isQuadratic = ( link2Nb->second < 0 );
               ++edgeVec[ isQuadratic ? SMDSEntity_Quad_Edge : SMDSEntity_Edge ];