Salome HOME
23126: [CEA 1562] Regression : Wrong nodes position using SetEnforcedVertex on a...
[modules/smesh.git] / src / SMESH / SMESH_MesherHelper.cxx
index 0556659dfe412cea035503dd439117b6c4d4529b..aec095ec3c4d00cfadc7dbad0dd626023c9befeb 100644 (file)
@@ -1,4 +1,4 @@
-// Copyright (C) 2007-2014  CEA/DEN, EDF R&D, OPEN CASCADE
+// Copyright (C) 2007-2015  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
@@ -32,6 +32,7 @@
 #include "SMDS_IteratorOnIterators.hxx"
 #include "SMDS_VolumeTool.hxx"
 #include "SMESH_Block.hxx"
+#include "SMESH_HypoFilter.hxx"
 #include "SMESH_MeshAlgos.hxx"
 #include "SMESH_ProxyMesh.hxx"
 #include "SMESH_subMesh.hxx"
@@ -70,7 +71,7 @@ using namespace std;
 
 namespace {
 
-  gp_XYZ XYZ(const SMDS_MeshNode* n) { return gp_XYZ(n->X(), n->Y(), n->Z()); }
+  inline SMESH_TNodeXYZ XYZ(const SMDS_MeshNode* n) { return SMESH_TNodeXYZ(n); }
 
   enum { U_periodic = 1, V_periodic = 2 };
 }
@@ -190,12 +191,13 @@ bool SMESH_MesherHelper::IsQuadraticSubMesh(const TopoDS_Shape& aSh)
     }
   }
 
-  if ( nbOldLinks == myTLinkNodeMap.size() )
+  // if ( nbOldLinks == myTLinkNodeMap.size() ) -- 0023068
+  if ( myTLinkNodeMap.empty() )
     myCreateQuadratic = false;
 
-  if(!myCreateQuadratic) {
+  if ( !myCreateQuadratic )
     myTLinkNodeMap.clear();
-  }
+
   SetSubShape( aSh );
 
   return myCreateQuadratic;
@@ -302,10 +304,20 @@ void SMESH_MesherHelper::SetSubShape(const TopoDS_Shape& aSh)
             isSeam = ( Abs( uv1.Coord(2) - myPar1[1] ) < Precision::PConfusion() ||
                        Abs( uv1.Coord(2) - myPar2[1] ) < Precision::PConfusion() );
           }
+          if ( isSeam ) // vertices are on period boundary, check a middle point (23032)
+          {
+            double f,l, r = 0.2345;
+            Handle(Geom2d_Curve) C2d = BRep_Tool::CurveOnSurface( edge, face, f, l );
+            uv2 = C2d->Value( f * r + l * ( 1.-r ));
+            if ( du < Precision::PConfusion() )
+              isSeam = ( Abs( uv1.Coord(1) - uv2.Coord(1) ) < Precision::PConfusion() );
+            else
+              isSeam = ( Abs( uv1.Coord(2) - uv2.Coord(2) ) < Precision::PConfusion() );
+          }
         }
         if ( isSeam )
         {
-          // store seam shape indices, negative if shape encounters twice
+          // store seam shape indices, negative if shape encounters twice ('real seam')
           mySeamShapeIds.insert( IsSeamShape( edgeID ) ? -edgeID : edgeID );
           for ( TopExp_Explorer v( edge, TopAbs_VERTEX ); v.More(); v.Next() ) {
             int vertexID = meshDS->ShapeToIndex( v.Current() );
@@ -803,18 +815,20 @@ namespace
 }
 
 //=======================================================================
-//function : applyIn2D
+//function : ApplyIn2D
 //purpose  : Perform given operation on two 2d points in parameric space of given surface.
 //           It takes into account period of the surface. Use gp_XY_FunPtr macro
 //           to easily define pointer to function of gp_XY class.
 //=======================================================================
 
-gp_XY SMESH_MesherHelper::applyIn2D(const Handle(Geom_Surface)& surface,
-                                    const gp_XY&                uv1,
-                                    const gp_XY&                uv2,
-                                    xyFunPtr                    fun,
-                                    const bool                  resultInPeriod)
+gp_XY SMESH_MesherHelper::ApplyIn2D(Handle(Geom_Surface) surface,
+                                    const gp_XY&         uv1,
+                                    const gp_XY&         uv2,
+                                    xyFunPtr             fun,
+                                    const bool           resultInPeriod)
 {
+  if ( surface->IsKind(STANDARD_TYPE(Geom_RectangularTrimmedSurface )))
+    surface = Handle(Geom_RectangularTrimmedSurface)::DownCast( surface )->BasisSurface();
   Standard_Boolean isUPeriodic = surface.IsNull() ? false : surface->IsUPeriodic();
   Standard_Boolean isVPeriodic = surface.IsNull() ? false : surface->IsVPeriodic();
   if ( !isUPeriodic && !isVPeriodic )
@@ -842,6 +856,31 @@ gp_XY SMESH_MesherHelper::applyIn2D(const Handle(Geom_Surface)& surface,
 
   return res;
 }
+
+//=======================================================================
+//function : AdjustByPeriod
+//purpose  : Move node positions on a FACE within surface period
+//=======================================================================
+
+void SMESH_MesherHelper::AdjustByPeriod( const TopoDS_Face& face, gp_XY uv[], const int nbUV )
+{
+  SMESH_MesherHelper h( *myMesh ), *ph = face.IsSame( myShape ) ? this : &h;
+  ph->SetSubShape( face );
+
+  for ( int iCoo = U_periodic; iCoo <= V_periodic; ++iCoo )
+    if ( ph->GetPeriodicIndex() & iCoo )
+    {
+      const double period = ( ph->myPar2[iCoo-1] - ph->myPar1[iCoo-1] );
+      const double xRef = uv[0].Coord( iCoo );
+      for ( int i = 1; i < nbUV; ++i )
+      {
+        double x = uv[i].Coord( iCoo );
+        double dx = ShapeAnalysis::AdjustByPeriod( x, xRef, period );
+        uv[i].SetCoord( iCoo, x + dx );
+      }
+    }
+}
+
 //=======================================================================
 //function : GetMiddleUV
 //purpose  : Return middle UV taking in account surface period
@@ -852,13 +891,13 @@ gp_XY SMESH_MesherHelper::GetMiddleUV(const Handle(Geom_Surface)& surface,
                                       const gp_XY&                p2)
 {
   // NOTE:
-  // the proper place of getting basic surface seems to be in applyIn2D()
+  // the proper place of getting basic surface seems to be in ApplyIn2D()
   // but we put it here to decrease a risk of regressions just before releasing a version
-  Handle(Geom_Surface) surf = surface;
-  while ( !surf.IsNull() && surf->IsKind(STANDARD_TYPE(Geom_RectangularTrimmedSurface )))
-    surf = Handle(Geom_RectangularTrimmedSurface)::DownCast( surf )->BasisSurface();
+  // Handle(Geom_Surface) surf = surface;
+  // while ( !surf.IsNull() && surf->IsKind(STANDARD_TYPE(Geom_RectangularTrimmedSurface )))
+  //   surf = Handle(Geom_RectangularTrimmedSurface)::DownCast( surf )->BasisSurface();
 
-  return applyIn2D( surf, p1, p2, & AverageUV );
+  return ApplyIn2D( surface, p1, p2, & AverageUV );
 }
 
 //=======================================================================
@@ -1256,14 +1295,24 @@ const SMDS_MeshNode* SMESH_MesherHelper::GetCentralNode(const SMDS_MeshNode* n1,
   if ( solidID < 1 && !faceId2nbNodes.empty() ) // SOLID not found
   {
     // find ID of the FACE the four corner nodes belong to
-    itMapWithIdFace = faceId2nbNodes.begin();
-    for ( ; itMapWithIdFace != faceId2nbNodes.end(); ++itMapWithIdFace)
+    itMapWithIdFace = faceId2nbNodes.find( myShapeID ); // IPAL52698
+    if ( itMapWithIdFace != faceId2nbNodes.end() &&
+         itMapWithIdFace->second == 4 )
+    {
+      shapeType = TopAbs_FACE;
+      faceID = myShapeID;
+    }
+    else
     {
-      if ( itMapWithIdFace->second == 4 ) 
+      itMapWithIdFace = faceId2nbNodes.begin();
+      for ( ; itMapWithIdFace != faceId2nbNodes.end(); ++itMapWithIdFace)
       {
-        shapeType = TopAbs_FACE;
-        faceID = (*itMapWithIdFace).first;
-        break;
+        if ( itMapWithIdFace->second == 4 )
+        {
+          shapeType = TopAbs_FACE;
+          faceID = (*itMapWithIdFace).first;
+          break;
+        }
       }
     }
   }
@@ -1281,11 +1330,20 @@ const SMDS_MeshNode* SMESH_MesherHelper::GetCentralNode(const SMDS_MeshNode* n1,
   bool toCheck = true;
   if ( !F.IsNull() && !force3d )
   {
-    uvAvg = calcTFI (0.5, 0.5,
-                     GetNodeUV(F,n1,n3,&toCheck), GetNodeUV(F,n2,n4,&toCheck),
-                     GetNodeUV(F,n3,n1,&toCheck), GetNodeUV(F,n4,n2,&toCheck), 
-                     GetNodeUV(F,n12,n3), GetNodeUV(F,n23,n4),
-                     GetNodeUV(F,n34,n2), GetNodeUV(F,n41,n2));
+    gp_XY uv[8] = {
+      GetNodeUV( F,n1,  n3, &toCheck ),
+      GetNodeUV( F,n2,  n4, &toCheck ),
+      GetNodeUV( F,n3,  n1, &toCheck ),
+      GetNodeUV( F,n4,  n2, &toCheck ),
+      GetNodeUV( F,n12, n3 ),
+      GetNodeUV( F,n23, n4 ),
+      GetNodeUV( F,n34, n2 ),
+      GetNodeUV( F,n41, n2 )
+    };
+    AdjustByPeriod( F, uv, 8 ); // put uv[] within a period (IPAL52698)
+
+    uvAvg = calcTFI (0.5, 0.5, uv[0],uv[1],uv[2],uv[3], uv[4],uv[5],uv[6],uv[7] );
+
     TopLoc_Location loc;
     Handle( Geom_Surface ) S = BRep_Tool::Surface( F, loc );
     P = S->Value( uvAvg.X(), uvAvg.Y() ).Transformed( loc );
@@ -1297,7 +1355,7 @@ const SMDS_MeshNode* SMESH_MesherHelper::GetCentralNode(const SMDS_MeshNode* n1,
   {
     P = calcTFI (0.5, 0.5,
                  SMESH_TNodeXYZ(n1),  SMESH_TNodeXYZ(n2),
-                 SMESH_TNodeXYZ(n3),  SMESH_TNodeXYZ(n4), 
+                 SMESH_TNodeXYZ(n3),  SMESH_TNodeXYZ(n4),
                  SMESH_TNodeXYZ(n12), SMESH_TNodeXYZ(n23),
                  SMESH_TNodeXYZ(n34), SMESH_TNodeXYZ(n41));
     centralNode = meshDS->AddNode( P.X(), P.Y(), P.Z() );
@@ -1401,14 +1459,24 @@ const SMDS_MeshNode* SMESH_MesherHelper::GetCentralNode(const SMDS_MeshNode* n1,
   if ( solidID < 1 && !faceId2nbNodes.empty() ) // SOLID not found
   {
     // find ID of the FACE the four corner nodes belong to
-    itMapWithIdFace = faceId2nbNodes.begin();
-    for ( ; itMapWithIdFace != faceId2nbNodes.end(); ++itMapWithIdFace)
+    itMapWithIdFace = faceId2nbNodes.find( myShapeID ); // IPAL52698
+    if ( itMapWithIdFace != faceId2nbNodes.end() &&
+         itMapWithIdFace->second == 4 )
+    {
+      shapeType = TopAbs_FACE;
+      faceID = myShapeID;
+    }
+    else
     {
-      if ( itMapWithIdFace->second == 3 ) 
+      itMapWithIdFace = faceId2nbNodes.begin();
+      for ( ; itMapWithIdFace != faceId2nbNodes.end(); ++itMapWithIdFace)
       {
-        shapeType = TopAbs_FACE;
-        faceID = (*itMapWithIdFace).first;
-        break;
+        if ( itMapWithIdFace->second == 3 )
+        {
+          shapeType = TopAbs_FACE;
+          faceID = (*itMapWithIdFace).first;
+          break;
+        }
       }
     }
   }
@@ -1420,13 +1488,18 @@ const SMDS_MeshNode* SMESH_MesherHelper::GetCentralNode(const SMDS_MeshNode* n1,
   {
     F = TopoDS::Face( meshDS->IndexToShape( faceID ));
     bool checkOK = true, badTria = false;
-    gp_XY uv1  = GetNodeUV( F, n1, n23, &checkOK );
-    gp_XY uv2  = GetNodeUV( F, n2, n31, &checkOK );
-    gp_XY uv3  = GetNodeUV( F, n3, n12, &checkOK );
-    gp_XY uv12 = GetNodeUV( F, n12, n3, &checkOK );
-    gp_XY uv23 = GetNodeUV( F, n23, n1, &checkOK );
-    gp_XY uv31 = GetNodeUV( F, n31, n2, &checkOK );
-    uvAvg = GetCenterUV( uv1,uv2,uv3, uv12,uv23,uv31, &badTria );
+    gp_XY uv[6] = {
+      GetNodeUV( F, n1, n23, &checkOK ),
+      GetNodeUV( F, n2, n31, &checkOK ),
+      GetNodeUV( F, n3, n12, &checkOK ),
+      GetNodeUV( F, n12, n3, &checkOK ),
+      GetNodeUV( F, n23, n1, &checkOK ),
+      GetNodeUV( F, n31, n2, &checkOK )
+    };
+    AdjustByPeriod( F, uv, 6 ); // put uv[] within a period (IPAL52698)
+
+    uvAvg = GetCenterUV( uv[0],uv[1],uv[2], uv[3],uv[4],uv[5], &badTria );
+
     if ( badTria || !checkOK )
       force3d = true;
   }
@@ -1924,26 +1997,28 @@ SMDS_MeshFace* SMESH_MesherHelper::AddPolygonalFace (const vector<const SMDS_Mes
   SMESHDS_Mesh * meshDS = GetMeshDS();
   SMDS_MeshFace* elem = 0;
 
-  if(!myCreateQuadratic) {
+  if(!myCreateQuadratic)
+  {
     if(id)
       elem = meshDS->AddPolygonalFaceWithID(nodes, id);
     else
       elem = meshDS->AddPolygonalFace(nodes);
   }
-  else {
-    vector<const SMDS_MeshNode*> newNodes;
+  else
+  {
+    vector<const SMDS_MeshNode*> newNodes( nodes.size() * 2 );
+    newNodes = nodes;
     for ( int i = 0; i < nodes.size(); ++i )
     {
       const SMDS_MeshNode* n1 = nodes[i];
       const SMDS_MeshNode* n2 = nodes[(i+1)%nodes.size()];
       const SMDS_MeshNode* n12 = GetMediumNode( n1, n2, force3d, TopAbs_FACE );
-      newNodes.push_back( n1 );
       newNodes.push_back( n12 );
     }
     if(id)
-      elem = meshDS->AddPolygonalFaceWithID(newNodes, id);
+      elem = meshDS->AddQuadPolygonalFaceWithID(newNodes, id);
     else
-      elem = meshDS->AddPolygonalFace(newNodes);
+      elem = meshDS->AddQuadPolygonalFace(newNodes);
   }
   if ( mySetElemOnShape && myShapeID > 0 )
     meshDS->SetMeshElementOnShape( elem, myShapeID );
@@ -2125,7 +2200,7 @@ SMDS_MeshVolume* SMESH_MesherHelper::AddVolume(const SMDS_MeshNode* n1,
     const SMDS_MeshNode* n26 = GetMediumNode( n2, n6, force3d, TopAbs_SOLID );
     const SMDS_MeshNode* n37 = GetMediumNode( n3, n7, force3d, TopAbs_SOLID );
     const SMDS_MeshNode* n48 = GetMediumNode( n4, n8, force3d, TopAbs_SOLID );
-    if(myCreateBiQuadratic)
+    if ( myCreateBiQuadratic )
     {
       const SMDS_MeshNode* n1234 = GetCentralNode( n1,n2,n3,n4,n12,n23,n34,n41,force3d );
       const SMDS_MeshNode* n1256 = GetCentralNode( n1,n2,n5,n6,n12,n26,n56,n15,force3d );
@@ -2160,9 +2235,9 @@ SMDS_MeshVolume* SMESH_MesherHelper::AddVolume(const SMDS_MeshNode* n1,
 
       pointsOnShapes[ SMESH_Block::ID_Fxy0 ] = SMESH_TNodeXYZ( n3478 );
       pointsOnShapes[ SMESH_Block::ID_Fxy1 ] = SMESH_TNodeXYZ( n1256 );
-      pointsOnShapes[ SMESH_Block::ID_Fx0z ] = SMESH_TNodeXYZ( n1458 );   
-      pointsOnShapes[ SMESH_Block::ID_Fx1z ] = SMESH_TNodeXYZ( n2367 );   
-      pointsOnShapes[ SMESH_Block::ID_F0yz ] = SMESH_TNodeXYZ( n1234 );    
+      pointsOnShapes[ SMESH_Block::ID_Fx0z ] = SMESH_TNodeXYZ( n1458 );
+      pointsOnShapes[ SMESH_Block::ID_Fx1z ] = SMESH_TNodeXYZ( n2367 );
+      pointsOnShapes[ SMESH_Block::ID_F0yz ] = SMESH_TNodeXYZ( n1234 );
       pointsOnShapes[ SMESH_Block::ID_F1yz ] = SMESH_TNodeXYZ( n5678 );
 
       gp_XYZ centerCube(0.5, 0.5, 0.5);
@@ -2172,27 +2247,27 @@ SMDS_MeshVolume* SMESH_MesherHelper::AddVolume(const SMDS_MeshNode* n1,
         meshDS->AddNode( nCenterElem.X(), nCenterElem.Y(), nCenterElem.Z() );
       meshDS->SetNodeInVolume( nCenter, myShapeID );
 
-     if(id)
+      if(id)
         elem = meshDS->AddVolumeWithID(n1, n2, n3, n4, n5, n6, n7, n8,
-                                      n12, n23, n34, n41, n56, n67,
-                                      n78, n85, n15, n26, n37, n48,
-                                      n1234, n1256, n2367, n3478, n1458, n5678, nCenter, id);
+                                       n12, n23, n34, n41, n56, n67,
+                                       n78, n85, n15, n26, n37, n48,
+                                       n1234, n1256, n2367, n3478, n1458, n5678, nCenter, id);
       else
         elem = meshDS->AddVolume(n1, n2, n3, n4, n5, n6, n7, n8,
-                                n12, n23, n34, n41, n56, n67,
-                                n78, n85, n15, n26, n37, n48,
-                                n1234, n1256, n2367, n3478, n1458, n5678, nCenter);
+                                 n12, n23, n34, n41, n56, n67,
+                                 n78, n85, n15, n26, n37, n48,
+                                 n1234, n1256, n2367, n3478, n1458, n5678, nCenter);
     }
     else
     {
       if(id)
         elem = meshDS->AddVolumeWithID(n1, n2, n3, n4, n5, n6, n7, n8,
-                                      n12, n23, n34, n41, n56, n67,
-                                      n78, n85, n15, n26, n37, n48, id);
+                                       n12, n23, n34, n41, n56, n67,
+                                       n78, n85, n15, n26, n37, n48, id);
       else
         elem = meshDS->AddVolume(n1, n2, n3, n4, n5, n6, n7, n8,
-                                n12, n23, n34, n41, n56, n67,
-                                n78, n85, n15, n26, n37, n48);
+                                 n12, n23, n34, n41, n56, n67,
+                                 n78, n85, n15, n26, n37, n48);
     }
   }
   if ( mySetElemOnShape && myShapeID > 0 )
@@ -2571,7 +2646,7 @@ bool SMESH_MesherHelper::IsStructured( SMESH_subMesh* faceSM )
   faceAnalyser.SetSubShape( faceSM->GetSubShape() );
 
   // rotate edges to get the first node being at corner
-  // (in principle it's not necessary but so far none SALOME algo can make
+  // (in principle it's not necessary because so far none SALOME algo can make
   //  such a structured mesh that all corner nodes are not on VERTEXes)
   bool isCorner     = false;
   int nbRemainEdges = nbEdgesInWires.front();
@@ -2733,38 +2808,80 @@ bool SMESH_MesherHelper::IsReversedSubMesh (const TopoDS_Face& theFace)
   if ( !aSubMeshDSFace )
     return isReversed;
 
-  // find an element with a good normal
-  gp_Vec Ne;
-  bool normalOK = false;
-  gp_XY uv;
+  // find an element on a bounday of theFace
   SMDS_ElemIteratorPtr iteratorElem = aSubMeshDSFace->GetElements();
-  while ( !normalOK && iteratorElem->more() ) // loop on elements on theFace
+  const SMDS_MeshNode* nn[2];
+  while ( iteratorElem->more() ) // loop on elements on theFace
   {
     const SMDS_MeshElement* elem = iteratorElem->next();
-    if ( elem && elem->NbCornerNodes() > 2 )
+    if ( ! elem ) continue;
+
+    // look for 2 nodes on EDGE
+    int nbNodes = elem->NbCornerNodes();
+    nn[0] = elem->GetNode( nbNodes-1 );
+    for ( int iN = 0; iN < nbNodes; ++iN )
     {
-      SMESH_TNodeXYZ nPnt[3];
-      SMDS_ElemIteratorPtr nodesIt = elem->nodesIterator();
-      int iNodeOnFace = 0, iPosDim = SMDS_TOP_VERTEX;
-      for ( int iN = 0; nodesIt->more() && iN < 3; ++iN) // loop on nodes
+      nn[1] = elem->GetNode( iN );
+      if ( nn[0]->GetPosition()->GetDim() < 2 &&
+           nn[1]->GetPosition()->GetDim() < 2 )
       {
-        nPnt[ iN ] = nodesIt->next();
-        if ( nPnt[ iN ]._node->GetPosition()->GetTypeOfPosition() > iPosDim )
+        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 ))
         {
-          iNodeOnFace = iN;
-          iPosDim = nPnt[ iN ]._node->GetPosition()->GetTypeOfPosition();
+          // is E seam edge?
+          int nb = 0;
+          for ( TopExp_Explorer exp( theFace, TopAbs_EDGE ); exp.More(); exp.Next() )
+            if ( E.IsSame( exp.Current() )) {
+              ++nb;
+              E = exp.Current(); // to know orientation
+            }
+          if ( nb == 1 )
+          {
+            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 );
+            if ( ok )
+            {
+              isReversed = ( u0 > u1 );
+              if ( E.Orientation() == TopAbs_REVERSED )
+                isReversed = !isReversed;
+              return isReversed;
+            }
+          }
         }
       }
-      // compute normal
-      gp_Vec v01( nPnt[0], nPnt[1] ), v02( nPnt[0], nPnt[2] );
-      if ( v01.SquareMagnitude() > RealSmall() &&
-           v02.SquareMagnitude() > RealSmall() )
+      nn[0] = nn[1];
+    }
+  }
+
+  // find an element with a good normal
+  gp_Vec Ne;
+  bool normalOK = false;
+  gp_XY uv;
+  iteratorElem = aSubMeshDSFace->GetElements();
+  while ( !normalOK && iteratorElem->more() ) // loop on elements on theFace
+  {
+    const SMDS_MeshElement* elem = iteratorElem->next();
+    if ( ! SMESH_MeshAlgos::FaceNormal( elem, const_cast<gp_XYZ&>( Ne.XYZ() ), /*normalized=*/0 ))
+      continue;
+    normalOK = true;
+
+    // get UV of a node inside theFACE
+    SMDS_ElemIteratorPtr nodesIt = elem->nodesIterator();
+    const SMDS_MeshNode* nInFace = 0;
+    int iPosDim = SMDS_TOP_VERTEX;
+    while ( nodesIt->more() ) // loop on nodes
+    {
+      const SMDS_MeshNode* n = static_cast<const SMDS_MeshNode*>( nodesIt->next() );
+      if ( n->GetPosition()->GetTypeOfPosition() >= iPosDim )
       {
-        Ne = v01 ^ v02;
-        if (( normalOK = ( Ne.SquareMagnitude() > RealSmall() )))
-          uv = GetNodeUV( theFace, nPnt[iNodeOnFace]._node, 0, &normalOK );
+        nInFace = n;
+        iPosDim = n->GetPosition()->GetTypeOfPosition();
       }
     }
+    uv = GetNodeUV( theFace, nInFace, 0, &normalOK );
   }
   if ( !normalOK )
     return isReversed;
@@ -2775,11 +2892,8 @@ bool SMESH_MesherHelper::IsReversedSubMesh (const TopoDS_Face& theFace)
   // if ( surf.IsNull() || surf->Continuity() < GeomAbs_C1 )
   // some surfaces not detected as GeomAbs_C1 are nevertheless correct for meshing
   if ( surf.IsNull() || surf->Continuity() < GeomAbs_C0 )
-    {
-      if (!surf.IsNull())
-        MESSAGE("surf->Continuity() < GeomAbs_C1 " << (surf->Continuity() < GeomAbs_C1));
-      return isReversed;
-    }
+    return isReversed;
+
   gp_Vec d1u, d1v; gp_Pnt p;
   surf->D1( uv.X(), uv.Y(), p, d1u, d1v );
   gp_Vec Nf = (d1u ^ d1v).Transformed( loc );
@@ -3096,6 +3210,28 @@ TopAbs_ShapeEnum SMESH_MesherHelper::GetGroupType(const TopoDS_Shape& group,
   return TopAbs_SHAPE;
 }
 
+//================================================================================
+/*!
+ * \brief Returns a shape, to which a hypothesis used to mesh a given shape is assigned
+ *  \param [in] hyp - the hypothesis
+ *  \param [in] shape - the shape, for meshing which the \a hyp is used
+ *  \param [in] mesh - the mesh
+ *  \return TopoDS_Shape - the shape the \a hyp is assigned to
+ */
+//================================================================================
+
+TopoDS_Shape SMESH_MesherHelper::GetShapeOfHypothesis( const SMESHDS_Hypothesis * hyp,
+                                                       const TopoDS_Shape&        shape,
+                                                       SMESH_Mesh*                mesh)
+{
+  const SMESH_Hypothesis* h = static_cast<const SMESH_Hypothesis*>( hyp );
+  SMESH_HypoFilter hypFilter( SMESH_HypoFilter::Is( h ));
+
+  TopoDS_Shape shapeOfHyp;
+  mesh->GetHypothesis( shape, hypFilter, /*checkAncestors=*/true, &shapeOfHyp );
+  return shapeOfHyp;
+}
+
 //=======================================================================
 //function : IsQuadraticMesh
 //purpose  : Check mesh without geometry for: if all elements on this shape are quadratic,
@@ -3207,6 +3343,11 @@ TopoDS_Shape SMESH_MesherHelper::GetCommonAncestor(const TopoDS_Shape& shape1,
   TopoDS_Shape commonAnc;
   if ( !shape1.IsNull() && !shape2.IsNull() )
   {
+    if ( shape1.ShapeType() == ancestorType && IsSubShape( shape2, shape1 ))
+      return shape1;
+    if ( shape2.ShapeType() == ancestorType && IsSubShape( shape1, shape2 ))
+      return shape2;
+
     PShapeIteratorPtr ancIt = GetAncestors( shape1, mesh, ancestorType );
     while ( const TopoDS_Shape* anc = ancIt->next() )
       if ( IsSubShape( shape2, *anc ))
@@ -3247,30 +3388,32 @@ namespace { // Structures used by FixQuadraticElements()
     mutable vector<const QFace* > _faces;
     mutable gp_Vec                _nodeMove;
     mutable int                   _nbMoves;
+    mutable bool                  _is2dFixed; // is moved along surface or in 3D
 
     QLink(const SMDS_MeshNode* n1, const SMDS_MeshNode* n2, const SMDS_MeshNode* nm):
       SMESH_TLink( n1,n2 ), _mediumNode(nm), _nodeMove(0,0,0), _nbMoves(0) {
       _faces.reserve(4);
-      //if ( MediumPos() != SMDS_TOP_3DSPACE )
-        _nodeMove = MediumPnt() - MiddlePnt();
+      _nodeMove = MediumPnt() - MiddlePnt();
+      _is2dFixed = ( MediumPos() != SMDS_TOP_FACE );
     }
     void SetContinuesFaces() const;
     const QFace* GetContinuesFace( const QFace* face ) const;
-    bool OnBoundary() const;
+    bool   OnBoundary() const;
     gp_XYZ MiddlePnt() const { return ( XYZ( node1() ) + XYZ( node2() )) / 2.; }
     gp_XYZ MediumPnt() const { return XYZ( _mediumNode ); }
 
-    SMDS_TypeOfPosition MediumPos() const
+    SMDS_TypeOfPosition  MediumPos() const
     { return _mediumNode->GetPosition()->GetTypeOfPosition(); }
-    SMDS_TypeOfPosition EndPos(bool isSecond) const
+    SMDS_TypeOfPosition  EndPos(bool isSecond) const
     { return (isSecond ? node2() : node1())->GetPosition()->GetTypeOfPosition(); }
     const SMDS_MeshNode* EndPosNode(SMDS_TypeOfPosition pos) const
     { return EndPos(0) == pos ? node1() : EndPos(1) == pos ? node2() : 0; }
 
-    void Move(const gp_Vec& move, bool sum=false) const
-    { _nodeMove += move; _nbMoves += sum ? (_nbMoves==0) : 1; }
+    void Move(const gp_Vec& move, bool sum=false, bool is2dFixed=false) const
+    { _nodeMove += move; _nbMoves += sum ? (_nbMoves==0) : 1; _is2dFixed |= is2dFixed; }
     gp_XYZ Move() const { return _nodeMove.XYZ() / _nbMoves; }
     bool IsMoved() const { return (_nbMoves > 0 /*&& !IsStraight()*/); }
+    bool IsFixedOnSurface() const { return _is2dFixed; }
     bool IsStraight() const
     { return isStraightLink( (XYZ(node1())-XYZ(node2())).SquareModulus(),
                              _nodeMove.SquareMagnitude());
@@ -3749,7 +3892,7 @@ namespace { // Structures used by FixQuadraticElements()
       double r = thePrevLen / fullLen;
 
       gp_Vec move = linkNorm * refProj * ( 1 - r );
-      theLink->Move( move, true );
+      theLink->Move( move, /*sum=*/true );
 
       MSG(string(theStep,'.')<<" Move "<< theLink->_mediumNode->GetID()<<
           " by " << refProj * ( 1 - r ) << " following " <<
@@ -4740,7 +4883,7 @@ void SMESH_MesherHelper::FixQuadraticElements(SMESH_ComputeErrorPtr& compError,
         }
         else if ( error == ERR_TRI ) {  // chain contains continues triangles
           TSplitTriaResult res = splitTrianglesIntoChains( rawChain, chains, pos );
-          if ( res != _OK ) { // not quadrangles split into triangles
+          if ( res != _OK ) { // not 'quadrangles split into triangles' in chain
             fixTriaNearBoundary( rawChain, *this );
             break;
           }
@@ -4768,12 +4911,14 @@ void SMESH_MesherHelper::FixQuadraticElements(SMESH_ComputeErrorPtr& compError,
           // mesure chain length and compute link position along the chain
           double chainLen = 0;
           vector< double > linkPos;
+          TChain savedChain; // backup
           MSGBEG( "Link medium nodes: ");
           TChain::iterator link0 = chain.begin(), link1 = chain.begin(), link2;
           for ( ++link1; link1 != chain.end(); ++link1, ++link0 ) {
             MSGBEG( (*link0)->_mediumNode->GetID() << "-" <<(*link1)->_mediumNode->GetID()<<" ");
             double len = ((*link0)->MiddlePnt() - (*link1)->MiddlePnt()).Modulus();
             while ( len < numeric_limits<double>::min() ) { // remove degenerated link
+              if ( savedChain.empty() ) savedChain = chain;
               link1 = chain.erase( link1 );
               if ( link1 == chain.end() )
                 break;
@@ -4783,9 +4928,16 @@ void SMESH_MesherHelper::FixQuadraticElements(SMESH_ComputeErrorPtr& compError,
             linkPos.push_back( chainLen );
           }
           MSG("");
-          if ( linkPos.size() < 2 )
-            continue;
-
+          if ( linkPos.size() <= 2 && savedChain.size() > 2 ) {
+            //continue;
+            linkPos.clear();
+            chainLen = 0;
+            chain = savedChain;
+            for ( link1 = chain.begin(); link1 != chain.end(); ++link1 ) {
+              chainLen += 1;
+              linkPos.push_back( chainLen );
+            }
+          }
           gp_Vec move0 = chain.front()->_nodeMove;
           gp_Vec move1 = chain.back ()->_nodeMove;
 
@@ -4801,7 +4953,7 @@ void SMESH_MesherHelper::FixQuadraticElements(SMESH_ComputeErrorPtr& compError,
               face = TopoDS::Face( f );
               faceHlp.SetSubShape( face );
               Handle(Geom_Surface) surf = BRep_Tool::Surface(face,loc);
-              bool isStraight[2];
+              //bool isStraight[2]; // commented for issue 0023118
               for ( int is1 = 0; is1 < 2; ++is1 ) // move0 or move1
               {
                 TChainLink& link = is1 ? chain.back() : chain.front();
@@ -4810,17 +4962,17 @@ void SMESH_MesherHelper::FixQuadraticElements(SMESH_ComputeErrorPtr& compError,
                 gp_XY uv2  = faceHlp.GetNodeUV( face, link->node2(),     nodeOnFace, &checkUV );
                 gp_XY uv12 = faceHlp.GetMiddleUV( surf, uv1, uv2 );
                 // uvMove = uvm - uv12
-                gp_XY uvMove = applyIn2D(surf, uvm, uv12, gp_XY_Subtracted, /*inPeriod=*/false);
+                gp_XY uvMove = ApplyIn2D(surf, uvm, uv12, gp_XY_Subtracted, /*inPeriod=*/false);
                 ( is1 ? move1 : move0 ).SetCoord( uvMove.X(), uvMove.Y(), 0 );
                 if ( !is1 ) // correct nodeOnFace for move1 (issue 0020919)
                   nodeOnFace = (*(++chain.rbegin()))->_mediumNode;
-                isStraight[is1] = isStraightLink( (uv2-uv1).SquareModulus(),
-                                                  10 * uvMove.SquareModulus());
-              }
-              if ( isStraight[0] && isStraight[1] ) {
-                MSG("2D straight - ignore");
-                continue; // straight - no need to move nodes of internal links
+                // isStraight[is1] = isStraightLink( (uv2-uv1).SquareModulus(),
+                //                                   10 * uvMove.SquareModulus());
               }
+              // if ( isStraight[0] && isStraight[1] ) {
+              //   MSG("2D straight - ignore");
+              //   continue; // straight - no need to move nodes of internal links
+              // }
 
               // check if a chain is already fixed
               gp_XY uvm  = faceHlp.GetNodeUV( face, linkOnFace->_mediumNode, 0, &checkUV );
@@ -4864,15 +5016,20 @@ void SMESH_MesherHelper::FixQuadraticElements(SMESH_ComputeErrorPtr& compError,
               // transform to global
               gp_Vec x01( (*link0)->MiddlePnt(), (*link1)->MiddlePnt() );
               gp_Vec x12( (*link1)->MiddlePnt(), (*link2)->MiddlePnt() );
-              gp_Vec x = x01.Normalized() + x12.Normalized();
-              trsf.SetTransformation( gp_Ax3( gp::Origin(), link1->Normal(), x), gp_Ax3() );
+              try {
+                gp_Vec x = x01.Normalized() + x12.Normalized();
+                trsf.SetTransformation( gp_Ax3( gp::Origin(), link1->Normal(), x), gp_Ax3() );
+              } catch ( Standard_Failure ) {
+                trsf.Invert();
+              }
               move.Transform(trsf);
+              (*link1)->Move( move, /*sum=*/false, /*is2dFixed=*/false );
             }
             else {
               // compute 3D displacement by 2D one
               Handle(Geom_Surface) s = BRep_Tool::Surface(face,loc);
               gp_XY oldUV   = faceHlp.GetNodeUV( face, (*link1)->_mediumNode, 0, &checkUV );
-              gp_XY newUV   = applyIn2D( s, oldUV, gp_XY( move.X(),move.Y()), gp_XY_Added );
+              gp_XY newUV   = ApplyIn2D( s, oldUV, gp_XY( move.X(),move.Y()), gp_XY_Added );
               gp_Pnt newPnt = s->Value( newUV.X(), newUV.Y());
               move = gp_Vec( XYZ((*link1)->_mediumNode), newPnt.Transformed(loc) );
               if ( SMDS_FacePosition* nPos =
@@ -4891,8 +5048,8 @@ void SMESH_MesherHelper::FixQuadraticElements(SMESH_ComputeErrorPtr& compError,
                      "newUV: "<<newUV.X()<<", "<<newUV.Y()<<" \t");
               }
 #endif
+              (*link1)->Move( move, /*sum=*/false, /*is2dFixed=*/true );
             }
-            (*link1)->Move( move );
             MSG( "Move " << (*link1)->_mediumNode->GetID() << " following "
                  << chain.front()->_mediumNode->GetID() <<"-"
                  << chain.back ()->_mediumNode->GetID() <<
@@ -4911,11 +5068,28 @@ void SMESH_MesherHelper::FixQuadraticElements(SMESH_ComputeErrorPtr& compError,
   const bool toFixCentralNodes = ( myMesh->NbBiQuadQuadrangles() +
                                    myMesh->NbBiQuadTriangles() +
                                    myMesh->NbTriQuadraticHexas() );
+  double distXYZ[4];
+  faceHlp.ToFixNodeParameters( true );
 
   for ( pLink = links.begin(); pLink != links.end(); ++pLink ) {
     if ( pLink->IsMoved() )
     {
       gp_Pnt p = pLink->MiddlePnt() + pLink->Move();
+
+      // put on surface nodes on FACE but moved in 3D (23050)
+      if ( !pLink->IsFixedOnSurface() )
+      {
+        faceHlp.SetSubShape( pLink->_mediumNode->getshapeId() );
+        if ( faceHlp.GetSubShape().ShapeType() == TopAbs_FACE )
+        {
+          const_cast<SMDS_MeshNode*>( pLink->_mediumNode )->setXYZ( p.X(), p.Y(), p.Z());
+          p.Coord( distXYZ[1], distXYZ[2], distXYZ[3] );
+          gp_XY uv( Precision::Infinite(), 0 );
+          if ( faceHlp.CheckNodeUV( TopoDS::Face( faceHlp.GetSubShape() ), pLink->_mediumNode,
+                                    uv, /*tol=*/pLink->Move().Modulus(), /*force=*/true, distXYZ ))
+            p.SetCoord( distXYZ[1], distXYZ[2], distXYZ[3] );
+        }
+      }
       GetMeshDS()->MoveNode( pLink->_mediumNode, p.X(), p.Y(), p.Z());
 
       // collect bi-quadratic elements
@@ -4964,6 +5138,7 @@ void SMESH_MesherHelper::FixQuadraticElements(SMESH_ComputeErrorPtr& compError,
         if ( i > 3 && nodes[i]->GetPosition()->GetTypeOfPosition() == SMDS_TOP_FACE )
           CheckNodeUV( F, nodes[i], uv[ i ], 2*tol, /*force=*/true );
       }
+      AdjustByPeriod( F, uv, 8 ); // put uv[] within a period (IPAL52698)
       // move the central node
       gp_XY uvCent = calcTFI (0.5, 0.5, uv[0],uv[1],uv[2],uv[3],uv[4],uv[5],uv[6],uv[7] );
       gp_Pnt p = surf->Value( uvCent.X(), uvCent.Y() ).Transformed( loc );
@@ -4989,7 +5164,7 @@ void SMESH_MesherHelper::FixQuadraticElements(SMESH_ComputeErrorPtr& compError,
       // nodes
       nodes.assign( tria->begin_nodes(), tria->end_nodes() );
       // UV
-      bool uvOK = true, badTria;
+      bool uvOK = true, badTria = false;
       for ( int i = 0; i < 6; ++i )
       {
         uv[ i ] = GetNodeUV( F, nodes[i], nodes[(i+1)%3], &uvOK );
@@ -4998,6 +5173,7 @@ void SMESH_MesherHelper::FixQuadraticElements(SMESH_ComputeErrorPtr& compError,
         if ( nodes[i]->GetPosition()->GetTypeOfPosition() == SMDS_TOP_FACE )
           CheckNodeUV( F, nodes[i], uv[ i ], 2*tol, /*force=*/true );
       }
+
       // move the central node
       gp_Pnt p;
       if ( !uvOK || badTria )
@@ -5008,6 +5184,7 @@ void SMESH_MesherHelper::FixQuadraticElements(SMESH_ComputeErrorPtr& compError,
       }
       else
       {
+        AdjustByPeriod( F, uv, 6 ); // put uv[] within a period (IPAL52698)
         gp_XY uvCent = GetCenterUV( uv[0], uv[1], uv[2], uv[3], uv[4], uv[5], &badTria );
         p = surf->Value( uvCent.X(), uvCent.Y() ).Transformed( loc );
       }