Salome HOME
Update copyright
[modules/smesh.git] / src / SMESH / SMESH_MesherHelper.cxx
index 12232f0b3c80c3a946ad073f59f7c689fb8325d1..5d26bb47295a3363ba87d783c15ff36ac0545add 100644 (file)
@@ -1,23 +1,23 @@
-//  Copyright (C) 2007-2010  CEA/DEN, EDF R&D, OPEN CASCADE
+// Copyright (C) 2007-2011  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
 //
 
 // File:      SMESH_MesherHelper.cxx
@@ -40,6 +40,7 @@
 #include <GeomAPI_ProjectPointOnCurve.hxx>
 #include <GeomAPI_ProjectPointOnSurf.hxx>
 #include <Geom_Curve.hxx>
+//#include <Geom_RectangularTrimmedSurface.hxx>
 #include <Geom_Surface.hxx>
 #include <ShapeAnalysis.hxx>
 #include <TopExp.hxx>
@@ -123,35 +124,50 @@ bool SMESH_MesherHelper::IsQuadraticSubMesh(const TopoDS_Shape& aSh)
 
   int nbOldLinks = myTLinkNodeMap.size();
 
-  TopExp_Explorer exp( aSh, subType );
-  for (; exp.More() && myCreateQuadratic; exp.Next()) {
-    if ( SMESHDS_SubMesh * subMesh = meshDS->MeshElements( exp.Current() )) {
-      if ( SMDS_ElemIteratorPtr it = subMesh->GetElements() ) {
-        while(it->more()) {
-          const SMDS_MeshElement* e = it->next();
-          if ( e->GetType() != elemType || !e->IsQuadratic() ) {
-            myCreateQuadratic = false;
-            break;
-          }
-          else {
-            // fill TLinkNodeMap
-            switch ( e->NbNodes() ) {
-            case 3:
-              AddTLinkNode(e->GetNode(0),e->GetNode(1),e->GetNode(2)); break;
-            case 6:
-              AddTLinkNode(e->GetNode(0),e->GetNode(1),e->GetNode(3));
-              AddTLinkNode(e->GetNode(1),e->GetNode(2),e->GetNode(4));
-              AddTLinkNode(e->GetNode(2),e->GetNode(0),e->GetNode(5)); break;
-            case 8:
-              AddTLinkNode(e->GetNode(0),e->GetNode(1),e->GetNode(4));
-              AddTLinkNode(e->GetNode(1),e->GetNode(2),e->GetNode(5));
-              AddTLinkNode(e->GetNode(2),e->GetNode(3),e->GetNode(6));
-              AddTLinkNode(e->GetNode(3),e->GetNode(0),e->GetNode(7));
-              break;
-            default:
+  if ( !myMesh->HasShapeToMesh() )
+  {
+    if (( myCreateQuadratic = myMesh->NbFaces( ORDER_QUADRATIC )))
+    {
+      SMDS_FaceIteratorPtr fIt = meshDS->facesIterator();
+      while ( fIt->more() )
+        AddTLinks( static_cast< const SMDS_MeshFace* >( fIt->next() ));
+    }
+  }
+  else
+  {
+    TopExp_Explorer exp( aSh, subType );
+    TopTools_MapOfShape checkedSubShapes;
+    for (; exp.More() && myCreateQuadratic; exp.Next()) {
+      if ( !checkedSubShapes.Add( exp.Current() ))
+        continue; // needed if aSh is compound of solids
+      if ( SMESHDS_SubMesh * subMesh = meshDS->MeshElements( exp.Current() )) {
+        if ( SMDS_ElemIteratorPtr it = subMesh->GetElements() ) {
+          while(it->more()) {
+            const SMDS_MeshElement* e = it->next();
+            if ( e->GetType() != elemType || !e->IsQuadratic() ) {
               myCreateQuadratic = false;
               break;
             }
+            else {
+              // fill TLinkNodeMap
+              switch ( e->NbNodes() ) {
+              case 3:
+                AddTLinkNode(e->GetNode(0),e->GetNode(1),e->GetNode(2)); break;
+              case 6:
+                AddTLinkNode(e->GetNode(0),e->GetNode(1),e->GetNode(3));
+                AddTLinkNode(e->GetNode(1),e->GetNode(2),e->GetNode(4));
+                AddTLinkNode(e->GetNode(2),e->GetNode(0),e->GetNode(5)); break;
+              case 8:
+                AddTLinkNode(e->GetNode(0),e->GetNode(1),e->GetNode(4));
+                AddTLinkNode(e->GetNode(1),e->GetNode(2),e->GetNode(5));
+                AddTLinkNode(e->GetNode(2),e->GetNode(3),e->GetNode(6));
+                AddTLinkNode(e->GetNode(3),e->GetNode(0),e->GetNode(7));
+                break;
+              default:
+                myCreateQuadratic = false;
+                break;
+              }
+            }
           }
         }
       }
@@ -210,9 +226,15 @@ void SMESH_MesherHelper::SetSubShape(const TopoDS_Shape& aSh)
   for ( TopExp_Explorer eF( aSh, TopAbs_FACE ); eF.More(); eF.Next() )
   {
     const TopoDS_Face& face = TopoDS::Face( eF.Current() );
-    BRepAdaptor_Surface surface( face );
-    if ( surface.IsUPeriodic() || surface.IsVPeriodic() )
+    TopLoc_Location loc;
+    Handle(Geom_Surface) surface = BRep_Tool::Surface( face, loc );
+
+    if ( surface->IsUPeriodic() || surface->IsVPeriodic() )
     {
+      //while ( surface->IsKind(STANDARD_TYPE(Geom_RectangularTrimmedSurface )))
+      //surface = Handle(Geom_RectangularTrimmedSurface)::DownCast( surface )->BasisSurface();
+      GeomAdaptor_Surface surf( surface );
+
       for (TopExp_Explorer exp( face, TopAbs_EDGE ); exp.More(); exp.Next())
       {
         // look for a seam edge
@@ -224,13 +246,13 @@ void SMESH_MesherHelper::SetSubShape(const TopoDS_Shape& aSh)
           if ( Abs( uv1.Coord(1) - uv2.Coord(1) ) < Abs( uv1.Coord(2) - uv2.Coord(2) ))
           {
             myParIndex |= U_periodic;
-            myPar1[0] = surface.FirstUParameter();
-            myPar2[0] = surface.LastUParameter();
+            myPar1[0] = surf.FirstUParameter();
+            myPar2[0] = surf.LastUParameter();
           }
           else {
             myParIndex |= V_periodic;
-            myPar1[1] = surface.FirstVParameter();
-            myPar2[1] = surface.LastVParameter();
+            myPar1[1] = surf.FirstVParameter();
+            myPar2[1] = surf.LastVParameter();
           }
           // store seam shape indices, negative if shape encounters twice
           int edgeID = meshDS->ShapeToIndex( edge );
@@ -334,19 +356,21 @@ void SMESH_MesherHelper::AddTLinks(const SMDS_MeshEdge* edge)
  */
 //================================================================================
 
-void SMESH_MesherHelper::AddTLinks(const SMDS_MeshFace* face)
+void SMESH_MesherHelper::AddTLinks(const SMDS_MeshFace* f)
 {
-  if ( face->IsQuadratic() )
-  {
-    const int nbLinks = face->NbCornerNodes();
-    for ( int i = 0; i < nbLinks; ++i )
-    {
-      const SMDS_MeshNode* n1  = face->GetNode( i );
-      const SMDS_MeshNode* n2  = face->GetNode(( i + 1 ) % nbLinks );
-      const SMDS_MeshNode* n12 = face->GetNode( i + nbLinks );
-      AddTLinkNode( n1, n2, n12 );
+  if ( !f->IsPoly() )
+    switch ( f->NbNodes() ) {
+    case 6:
+      AddTLinkNode(f->GetNode(0),f->GetNode(1),f->GetNode(3));
+      AddTLinkNode(f->GetNode(1),f->GetNode(2),f->GetNode(4));
+      AddTLinkNode(f->GetNode(2),f->GetNode(0),f->GetNode(5)); break;
+    case 8:
+      AddTLinkNode(f->GetNode(0),f->GetNode(1),f->GetNode(4));
+      AddTLinkNode(f->GetNode(1),f->GetNode(2),f->GetNode(5));
+      AddTLinkNode(f->GetNode(2),f->GetNode(3),f->GetNode(6));
+      AddTLinkNode(f->GetNode(3),f->GetNode(0),f->GetNode(7));
+    default:;
     }
-  }
 }
 
 //================================================================================
@@ -473,7 +497,7 @@ gp_XY SMESH_MesherHelper::GetNodeUV(const TopoDS_Face&   F,
     if ( validU )
       uv = C2d->Value( u );
     else
-      uv.SetCoord(0.,0.);
+      uv.SetCoord( Precision::Infinite(),0.);
     if ( check || !validU )
       uvOK = CheckNodeUV( F, n, uv.ChangeCoord(), 10*MaxTolerance( F ),/*force=*/ !validU );
 
@@ -552,6 +576,10 @@ gp_XY SMESH_MesherHelper::GetNodeUV(const TopoDS_Face&   F,
         uv = GetUVOnSeam( uv, GetNodeUV( F, n2, 0 ));
     }
   }
+  else
+  {
+    uvOK = CheckNodeUV( F, n, uv.ChangeCoord(), 10*MaxTolerance( F ));
+  }
 
   if ( check )
     *check = uvOK;
@@ -572,18 +600,15 @@ bool SMESH_MesherHelper::CheckNodeUV(const TopoDS_Face&   F,
                                      double               distXYZ[4]) const
 {
   int shapeID = n->getshapeId();
-  if ( force || toCheckPosOnShape( shapeID ))
+  bool infinit = ( Precision::IsInfinite( uv.X() ) || Precision::IsInfinite( uv.Y() ));
+  if ( force || toCheckPosOnShape( shapeID ) || infinit )
   {
-    double toldis = tol;
-    double tolmin = 1.e-7*myMesh->GetMeshDS()->getMaxDim(); // nodes coordinates are stored in float format
-    if (toldis < tolmin) toldis = tolmin;
     // check that uv is correct
     TopLoc_Location loc;
     Handle(Geom_Surface) surface = BRep_Tool::Surface( F,loc );
     gp_Pnt nodePnt = XYZ( n ), surfPnt(0,0,0);
     double dist = 0;
     if ( !loc.IsIdentity() ) nodePnt.Transform( loc.Transformation().Inverted() );
-    bool infinit = ( Precision::IsInfinite( uv.X() ) || Precision::IsInfinite( uv.Y() ));
     if ( infinit ||
          (dist = nodePnt.Distance( surfPnt = surface->Value( uv.X(), uv.Y() ))) > tol )
     {
@@ -776,10 +801,6 @@ bool SMESH_MesherHelper::CheckNodeU(const TopoDS_Edge&   E,
   int shapeID = n->getshapeId();
   if ( force || toCheckPosOnShape( shapeID ))
   {
-    //double toldis = tol;
-    //double tolmin = 1.e-7*myMesh->GetMeshDS()->getMaxDim(); // nodes coordinates are stored in float format
-    //if (toldis < tolmin) toldis = tolmin;
-    // check that u is correct
     TopLoc_Location loc; double f,l;
     Handle(Geom_Curve) curve = BRep_Tool::Curve( E,loc,f,l );
     if ( curve.IsNull() ) // degenerated edge
@@ -801,7 +822,7 @@ bool SMESH_MesherHelper::CheckNodeU(const TopoDS_Edge&   E,
         distXYZ[0] = dist;
         distXYZ[1] = curvPnt.X(); distXYZ[2] = curvPnt.Y(); distXYZ[3]=curvPnt.Z();
       }
-      if ( dist > tol /*toldis*/ )
+      if ( dist > tol )
       {
         setPosOnShapeValidity( shapeID, false );
         // u incorrect, project the node to the curve
@@ -830,10 +851,10 @@ bool SMESH_MesherHelper::CheckNodeU(const TopoDS_Edge&   E,
           distXYZ[0] = dist;
           distXYZ[1] = curvPnt.X(); distXYZ[2] = curvPnt.Y(); distXYZ[3]=curvPnt.Z();
         }
-        if ( dist > tol /*toldis*/)
+        if ( dist > tol )
         {
           MESSAGE( "SMESH_MesherHelper::CheckNodeU(), invalid projection" );
-          MESSAGE("distance " << dist << " " << tol/*dis*/);
+          MESSAGE("distance " << dist << " " << tol );
           return false;
         }
         // store the fixed U on the edge
@@ -1742,6 +1763,28 @@ bool SMESH_MesherHelper::IsClosedEdge( const TopoDS_Edge& anEdge )
   return TopExp::FirstVertex( anEdge ).IsSame( TopExp::LastVertex( anEdge ));
 }
 
+//================================================================================
+/*!
+ * \brief Wrapper over TopExp::FirstVertex() and TopExp::LastVertex() fixing them
+ *  in the case of INTERNAL edge
+ */
+//================================================================================
+
+TopoDS_Vertex SMESH_MesherHelper::IthVertex( const bool  is2nd,
+                                             TopoDS_Edge anEdge,
+                                             const bool  CumOri )
+{
+  if ( anEdge.Orientation() >= TopAbs_INTERNAL )
+    anEdge.Orientation( TopAbs_FORWARD );
+
+  const TopAbs_Orientation tgtOri = is2nd ? TopAbs_REVERSED : TopAbs_FORWARD;
+  TopoDS_Iterator vIt( anEdge, CumOri );
+  while ( vIt.More() && vIt.Value().Orientation() != tgtOri )
+    vIt.Next();
+
+  return ( vIt.More() ? TopoDS::Vertex(vIt.Value()) : TopoDS_Vertex() );
+}
+
 //=======================================================================
 //function : IsQuadraticMesh
 //purpose  : Check mesh without geometry for: if all elements on this shape are quadratic,
@@ -2736,7 +2779,9 @@ void SMESH_MesherHelper::FixQuadraticElements(bool volumeOnly)
           faces.Add( f.Current() ); // in not meshed solid
       }
       else { // fix nodes in the solid and its faces
+#ifdef _DEBUG_
         MSG("FIX SOLID " << nbSolids-- << " #" << GetMeshDS()->ShapeToIndex(s.Current()));
+#endif
         SMESH_MesherHelper h(*myMesh);
         h.SetSubShape( s.Current() );
         h.FixQuadraticElements(false);