Salome HOME
IPAL52499: Prismatic mesh is not computed on a prismatic shape
[modules/smesh.git] / src / StdMeshers / StdMeshers_Prism_3D.cxx
index 6227072d9ab3c5396b771dff47f6b624ffba0731..5e253829ab1af52a83a2b8a25489a086df85af4a 100644 (file)
@@ -1,4 +1,4 @@
-// Copyright (C) 2007-2015  CEA/DEN, EDF R&D, OPEN CASCADE
+// Copyright (C) 2007-2016  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
@@ -51,6 +51,7 @@
 #include <Geom2d_Line.hxx>
 #include <GeomLib_IsPlanarSurface.hxx>
 #include <Geom_Curve.hxx>
+#include <TColStd_DataMapOfIntegerInteger.hxx>
 #include <TopExp.hxx>
 #include <TopExp_Explorer.hxx>
 #include <TopTools_ListIteratorOfListOfShape.hxx>
@@ -62,6 +63,7 @@
 #include <gp_Ax3.hxx>
 
 #include <limits>
+#include <numeric>
 
 using namespace std;
 
@@ -161,6 +163,12 @@ namespace {
     {
       return _src2tgtNodes;
     }
+    void SetEventListener( SMESH_subMesh* tgtSubMesh )
+    {
+      NSProjUtils::SetEventListener( tgtSubMesh,
+                                     _sourceHypo->GetSourceFace(),
+                                     _sourceHypo->GetSourceMesh() );
+    }
   };
   //=======================================================================
   /*!
@@ -339,7 +347,7 @@ namespace {
     // gravity center of a layer
     gp_XYZ O(0,0,0);
     int vertexCol = -1;
-    for ( int i = 0; i < columns.size(); ++i )
+    for ( size_t i = 0; i < columns.size(); ++i )
     {
       O += gpXYZ( (*columns[ i ])[ z ]);
       if ( vertexCol < 0 &&
@@ -351,7 +359,7 @@ namespace {
     // Z axis
     gp_Vec Z(0,0,0);
     int iPrev = columns.size()-1;
-    for ( int i = 0; i < columns.size(); ++i )
+    for ( size_t i = 0; i < columns.size(); ++i )
     {
       gp_Vec v1( O, gpXYZ( (*columns[ iPrev ])[ z ]));
       gp_Vec v2( O, gpXYZ( (*columns[ i ]    )[ z ]));
@@ -363,11 +371,11 @@ namespace {
     {
       O = gpXYZ( (*columns[ vertexCol ])[ z ]);
     }
-    if ( xColumn < 0 || xColumn >= columns.size() )
+    if ( xColumn < 0 || xColumn >= (int) columns.size() )
     {
       // select a column for X dir
       double maxDist = 0;
-      for ( int i = 0; i < columns.size(); ++i )
+      for ( size_t i = 0; i < columns.size(); ++i )
       {
         double dist = ( O - gpXYZ((*columns[ i ])[ z ])).SquareModulus();
         if ( dist > maxDist )
@@ -454,9 +462,9 @@ namespace {
     std::advance( edgeIt, nbEdges-1 );
     TopoDS_Edge   prevE = *edgeIt;
     // bool isPrevStraight = SMESH_Algo::IsStraight( prevE );
-    int           iPrev = nbEdges - 1;
+    // int           iPrev = nbEdges - 1;
 
-    int iUnite = -1; // the first of united EDGEs
+    // int iUnite = -1; // the first of united EDGEs
 
     // analyse angles between EDGEs
     int nbCorners = 0;
@@ -524,7 +532,7 @@ namespace {
   void pointsToPython(const std::vector<gp_XYZ>& p)
   {
 #ifdef _DEBUG_
-    for ( int i = SMESH_Block::ID_V000; i < p.size(); ++i )
+    for ( size_t i = SMESH_Block::ID_V000; i < p.size(); ++i )
     {
       cout << "mesh.AddNode( " << p[i].X() << ", "<< p[i].Y() << ", "<< p[i].Z() << ") # " << i <<" " ;
       SMESH_Block::DumpShapeID( i, cout ) << endl;
@@ -560,7 +568,9 @@ StdMeshers_Prism_3D::StdMeshers_Prism_3D(int hypId, int studyId, SMESH_Gen* gen)
 //================================================================================
 
 StdMeshers_Prism_3D::~StdMeshers_Prism_3D()
-{}
+{
+  pointsToPython( std::vector<gp_XYZ>() ); // avoid warning: pointsToPython defined but not used
+}
 
 //=======================================================================
 //function : CheckHypothesis
@@ -863,7 +873,7 @@ bool StdMeshers_Prism_3D::Compute(SMESH_Mesh& theMesh, const TopoDS_Shape& theSh
 
 
     // TODO. there are other ways to find out the source FACE:
-    // propagation, topological similarity, ect.
+    // propagation, topological similarity, etc...
 
     // simply try to mesh all not meshed SOLIDs
     if ( meshedFaces.empty() )
@@ -933,12 +943,12 @@ bool StdMeshers_Prism_3D::getWallFaces( Prism_3D::TPrismTopo & thePrism,
 
   list< TopoDS_Edge >::iterator edge = thePrism.myBottomEdges.begin();
   std::list< int >::iterator     nbE = thePrism.myNbEdgesInWires.begin();
+  std::list< int > nbQuadsPerWire;
   int iE = 0;
-  double f,l;
   while ( edge != thePrism.myBottomEdges.end() )
   {
     ++iE;
-    if ( BRep_Tool::Curve( *edge, f,l ).IsNull() )
+    if ( SMESH_Algo::isDegenerated( *edge ))
     {
       edge = thePrism.myBottomEdges.erase( edge );
       --iE;
@@ -946,12 +956,14 @@ bool StdMeshers_Prism_3D::getWallFaces( Prism_3D::TPrismTopo & thePrism,
     }
     else
     {
+      bool hasWallFace = false;
       TopTools_ListIteratorOfListOfShape faceIt( edgeToFaces.FindFromKey( *edge ));
       for ( ; faceIt.More(); faceIt.Next() )
       {
         const TopoDS_Face& face = TopoDS::Face( faceIt.Value() );
         if ( !thePrism.myBottom.IsSame( face ))
         {
+          hasWallFace = true;
           Prism_3D::TQuadList quadList( 1, quadAlgo->CheckNbEdges( *mesh, face ));
           if ( !quadList.back() )
             return toSM( error(TCom("Side face #") << shapeID( face )
@@ -970,12 +982,23 @@ bool StdMeshers_Prism_3D::getWallFaces( Prism_3D::TPrismTopo & thePrism,
           break;
         }
       }
-      ++edge;
+      if ( hasWallFace )
+      {
+        ++edge;
+      }
+      else // seam edge (IPAL53561)
+      {
+        edge = thePrism.myBottomEdges.erase( edge );
+        --iE;
+        --(*nbE);
+      }
     }
     if ( iE == *nbE )
     {
       iE = 0;
       ++nbE;
+      int nbQuadPrev = std::accumulate( nbQuadsPerWire.begin(), nbQuadsPerWire.end(), 0 );
+      nbQuadsPerWire.push_back( thePrism.myWallQuads.size() - nbQuadPrev );
     }
   }
 
@@ -987,12 +1010,14 @@ bool StdMeshers_Prism_3D::getWallFaces( Prism_3D::TPrismTopo & thePrism,
   // that is not so evident in case of several WIREs in the bottom FACE
   thePrism.myRightQuadIndex.clear();
   for ( size_t i = 0; i < thePrism.myWallQuads.size(); ++i )
-    thePrism.myRightQuadIndex.push_back( i+1 );
-  list< int >::iterator nbEinW = thePrism.myNbEdgesInWires.begin();
-  for ( int iLeft = 0; nbEinW != thePrism.myNbEdgesInWires.end(); ++nbEinW )
   {
-    thePrism.myRightQuadIndex[ iLeft + *nbEinW - 1 ] = iLeft; // 1st EDGE index of a current WIRE
-    iLeft += *nbEinW;
+    thePrism.myRightQuadIndex.push_back( i+1 ); // OK for all but the last EDGE of a WIRE
+  }
+  list< int >::iterator nbQinW = nbQuadsPerWire.begin();
+  for ( int iLeft = 0; nbQinW != nbQuadsPerWire.end(); ++nbQinW )
+  {
+    thePrism.myRightQuadIndex[ iLeft + *nbQinW - 1 ] = iLeft; // for the last EDGE of a WIRE
+    iLeft += *nbQinW;
   }
 
   while ( totalNbFaces - faceMap.Extent() > 2 )
@@ -1072,7 +1097,7 @@ bool StdMeshers_Prism_3D::getWallFaces( Prism_3D::TPrismTopo & thePrism,
   {
     // now only top and bottom FACEs are not in the faceMap
     faceMap.Add( thePrism.myBottom );
-    for ( TopExp_Explorer f( thePrism.myShape3D, TopAbs_FACE );f.More(); f.Next() )
+    for ( TopExp_Explorer f( thePrism.myShape3D, TopAbs_FACE ); f.More(); f.Next() )
       if ( !faceMap.Contains( f.Current() )) {
         thePrism.myTop = TopoDS::Face( f.Current() );
         break;
@@ -1110,7 +1135,7 @@ bool StdMeshers_Prism_3D::compute(const Prism_3D::TPrismTopo& thePrism)
       ( ! botSM->GetAlgo() ||
         ! _gen->Compute( *botSM->GetFather(), botSM->GetSubShape(), /*shapeOnly=*/true )))
     return error( COMPERR_BAD_INPUT_MESH,
-                  TCom( "No mesher defined to compute the face #")
+                  TCom( "No mesher defined to compute the base face #")
                   << shapeID( thePrism.myBottom ));
 
   // Make all side FACEs of thePrism meshed with quads
@@ -1149,7 +1174,7 @@ bool StdMeshers_Prism_3D::compute(const Prism_3D::TPrismTopo& thePrism)
 
   // Projections on the top and bottom faces are taken from nodes existing
   // on these faces; find correspondence between bottom and top nodes
-  myUseBlock = false;
+  myUseBlock = false; // is set to true if projection is done using "block approach"
   myBotToColumnMap.clear();
   if ( !assocOrProjBottom2Top( bottomToTopTrsf, thePrism ) ) // it also fills myBotToColumnMap
     return false;
@@ -1157,38 +1182,60 @@ bool StdMeshers_Prism_3D::compute(const Prism_3D::TPrismTopo& thePrism)
 
   // Create nodes inside the block
 
-  // use transformation (issue 0020680, IPAL0052499)
-  StdMeshers_Sweeper sweeper;
-  double tol;
-  bool allowHighBndError;
-
   if ( !myUseBlock )
   {
+    // use transformation (issue 0020680, IPAL0052499) or a "straight line" approach
+    StdMeshers_Sweeper sweeper;
+
     // load boundary nodes into sweeper
     bool dummy;
+    const SMDS_MeshNode* prevN0 = 0, *prevN1 = 0;
     list< TopoDS_Edge >::const_iterator edge = thePrism.myBottomEdges.begin();
     for ( ; edge != thePrism.myBottomEdges.end(); ++edge )
     {
       int edgeID = meshDS->ShapeToIndex( *edge );
       TParam2ColumnMap* u2col = const_cast<TParam2ColumnMap*>
         ( myBlock.GetParam2ColumnMap( edgeID, dummy ));
-      TParam2ColumnMap::iterator u2colIt = u2col->begin();
-      for ( ; u2colIt != u2col->end(); ++u2colIt )
+
+      TParam2ColumnMap::iterator u2colIt = u2col->begin(), u2colEnd = u2col->end();
+      const SMDS_MeshNode* n0 = u2colIt->second[0];
+      const SMDS_MeshNode* n1 = u2col->rbegin()->second[0];
+      if ( n0 == prevN0 || n0 == prevN1 ) ++u2colIt;
+      if ( n1 == prevN0 || n1 == prevN1 ) --u2colEnd;
+      prevN0 = n0; prevN1 = n1;
+
+      for ( ; u2colIt != u2colEnd; ++u2colIt )
         sweeper.myBndColumns.push_back( & u2colIt->second );
     }
-    // load node columns inside the bottom face
+    // load node columns inside the bottom FACE
     TNode2ColumnMap::iterator bot_column = myBotToColumnMap.begin();
+    sweeper.myIntColumns.reserve( myBotToColumnMap.size() );
     for ( ; bot_column != myBotToColumnMap.end(); ++bot_column )
       sweeper.myIntColumns.push_back( & bot_column->second );
 
-    tol = getSweepTolerance( thePrism );
-    allowHighBndError = !isSimpleBottom( thePrism );
-  }
+    myHelper->SetElementsOnShape( true );
 
-  if ( !myUseBlock && sweeper.ComputeNodes( *myHelper, tol, allowHighBndError ))
-  {
+    // If all "vertical" EDGEs are straight, then all nodes of an internal node column
+    // are located on a line connecting the top node and the bottom node.
+    bool isStrightColunm = allVerticalEdgesStraight( thePrism );
+    if ( !isStrightColunm )
+    {
+      double tol = getSweepTolerance( thePrism );
+      bool allowHighBndError = !isSimpleBottom( thePrism );
+      myUseBlock = !sweeper.ComputeNodes( *myHelper, tol, allowHighBndError );
+    }
+    else if ( sweeper.CheckSameZ() )
+    {
+      myUseBlock = !sweeper.ComputeNodesOnStraightSameZ( *myHelper );
+    }
+    else
+    {
+      myUseBlock = !sweeper.ComputeNodesOnStraight( *myHelper, thePrism.myBottom, thePrism.myTop );
+    }
+    myHelper->SetElementsOnShape( false );
   }
-  else // use block approach
+
+  if ( myUseBlock ) // use block approach
   {
     // loop on nodes inside the bottom face
     Prism_3D::TNode prevBNode;
@@ -1196,7 +1243,8 @@ bool StdMeshers_Prism_3D::compute(const Prism_3D::TPrismTopo& thePrism)
     for ( ; bot_column != myBotToColumnMap.end(); ++bot_column )
     {
       const Prism_3D::TNode& tBotNode = bot_column->first; // bottom TNode
-      if ( tBotNode.GetPositionType() != SMDS_TOP_FACE )
+      if ( tBotNode.GetPositionType() != SMDS_TOP_FACE &&
+           myBlock.HasNodeColumn( tBotNode.myNode ))
         continue; // node is not inside the FACE
 
       // column nodes; middle part of the column are zero pointers
@@ -1286,6 +1334,7 @@ bool StdMeshers_Prism_3D::compute(const Prism_3D::TPrismTopo& thePrism)
   if ( !smDS ) return toSM( error(COMPERR_BAD_INPUT_MESH, "Null submesh"));
 
   // loop on bottom mesh faces
+  vector< const TNodeColumn* > columns;
   SMDS_ElemIteratorPtr faceIt = smDS->GetElements();
   while ( faceIt->more() )
   {
@@ -1295,24 +1344,26 @@ bool StdMeshers_Prism_3D::compute(const Prism_3D::TPrismTopo& thePrism)
 
     // find node columns for each node
     int nbNodes = face->NbCornerNodes();
-    vector< const TNodeColumn* > columns( nbNodes );
+    columns.resize( nbNodes );
     for ( int i = 0; i < nbNodes; ++i )
     {
       const SMDS_MeshNode* n = face->GetNode( i );
-      if ( n->GetPosition()->GetTypeOfPosition() == SMDS_TOP_FACE ) {
+      columns[ i ] = NULL;
+
+      if ( n->GetPosition()->GetTypeOfPosition() != SMDS_TOP_FACE )
+        columns[ i ] = myBlock.GetNodeColumn( n );
+
+      if ( !columns[ i ] )
+      {
         TNode2ColumnMap::iterator bot_column = myBotToColumnMap.find( n );
         if ( bot_column == myBotToColumnMap.end() )
-          return toSM( error(TCom("No nodes found above node ") << n->GetID() ));
-        columns[ i ] = & bot_column->second;
-      }
-      else {
-        columns[ i ] = myBlock.GetNodeColumn( n );
-        if ( !columns[ i ] )
           return toSM( error(TCom("No side nodes found above node ") << n->GetID() ));
+        columns[ i ] = & bot_column->second;
       }
     }
     // create prisms
-    AddPrisms( columns, myHelper );
+    if ( !AddPrisms( columns, myHelper ))
+      return toSM( error("Different 'vertical' discretization"));
 
   } // loop on bottom mesh faces
 
@@ -1322,7 +1373,7 @@ bool StdMeshers_Prism_3D::compute(const Prism_3D::TPrismTopo& thePrism)
 
   // update state of sub-meshes (mostly in order to erase improper errors)
   SMESH_subMesh* sm = myHelper->GetMesh()->GetSubMesh( thePrism.myShape3D );
-  SMESH_subMeshIteratorPtr smIt = sm->getDependsOnIterator(/*includeSelf=*/false);
+  SMESH_subMeshIteratorPtr smIt = sm->getDependsOnIterator(/*includeSelf=*/true);
   while ( smIt->more() )
   {
     sm = smIt->next();
@@ -1575,7 +1626,7 @@ bool StdMeshers_Prism_3D::computeWalls(const Prism_3D::TPrismTopo& thePrism)
       if ( ! fSM->IsMeshComputed() )
       {
         // Top EDGEs must be projections from the bottom ones
-        // to compute stuctured quad mesh on wall FACEs
+        // to compute structured quad mesh on wall FACEs
         // ---------------------------------------------------
         const TopoDS_Edge& botE = (*quad)->side[ QUAD_BOTTOM_SIDE ].grid->Edge(0);
         const TopoDS_Edge& topE = (*quad)->side[ QUAD_TOP_SIDE    ].grid->Edge(0);
@@ -1810,17 +1861,20 @@ bool StdMeshers_Prism_3D::Evaluate(SMESH_Mesh&         theMesh,
  */
 //================================================================================
 
-void StdMeshers_Prism_3D::AddPrisms( vector<const TNodeColumn*> & columns,
+bool StdMeshers_Prism_3D::AddPrisms( vector<const TNodeColumn*> & columns,
                                      SMESH_MesherHelper*          helper)
 {
-  int nbNodes = columns.size();
-  int nbZ     = columns[0]->size();
-  if ( nbZ < 2 ) return;
+  size_t nbNodes = columns.size();
+  size_t nbZ     = columns[0]->size();
+  if ( nbZ < 2 ) return false;
+  for ( size_t i = 1; i < nbNodes; ++i )
+    if ( columns[i]->size() != nbZ )
+      return false;
 
   // find out orientation
   bool isForward = true;
   SMDS_VolumeTool vTool;
-  int z = 1;
+  size_t z = 1;
   switch ( nbNodes ) {
   case 3: {
     SMDS_VolumeOfNodes tmpPenta ( (*columns[0])[z-1], // bottom
@@ -1906,7 +1960,7 @@ void StdMeshers_Prism_3D::AddPrisms( vector<const TNodeColumn*> & columns,
     vector<const SMDS_MeshNode*> nodes( 2*nbNodes + 4*nbNodes);
     for ( z = 1; z < nbZ; ++z )
     {
-      for ( int i = 0; i < nbNodes; ++i ) {
+      for ( size_t i = 0; i < nbNodes; ++i ) {
         nodes[ i             ] = (*columns[ i ])[z+iBase1]; // bottom or top
         nodes[ 2*nbNodes-i-1 ] = (*columns[ i ])[z+iBase2]; // top or bottom
         // side
@@ -1920,6 +1974,8 @@ void StdMeshers_Prism_3D::AddPrisms( vector<const TNodeColumn*> & columns,
     }
 
   } // switch ( nbNodes )
+
+  return true;
 }
 
 //================================================================================
@@ -1975,7 +2031,7 @@ bool StdMeshers_Prism_3D::assocOrProjBottom2Top( const gp_Trsf & bottomToTopTrsf
     n2nMapPtr = & TProjction2dAlgo::instance( this )->GetNodesMap();
   }
 
-  if ( !n2nMapPtr || n2nMapPtr->size() < botSMDS->NbNodes() )
+  if ( !n2nMapPtr || (int) n2nMapPtr->size() < botSMDS->NbNodes() )
   {
     // associate top and bottom faces
     NSProjUtils::TShapeShapeMap shape2ShapeMap;
@@ -2053,7 +2109,8 @@ bool StdMeshers_Prism_3D::assocOrProjBottom2Top( const gp_Trsf & bottomToTopTrsf
   {
     const SMDS_MeshNode* botNode = bN_tN->first;
     const SMDS_MeshNode* topNode = bN_tN->second;
-    if ( botNode->GetPosition()->GetTypeOfPosition() != SMDS_TOP_FACE )
+    if ( botNode->GetPosition()->GetTypeOfPosition() != SMDS_TOP_FACE &&
+         myBlock.HasNodeColumn( botNode ))
       continue; // wall columns are contained in myBlock
     // create node column
     Prism_3D::TNode bN( botNode );
@@ -2250,29 +2307,39 @@ bool StdMeshers_Prism_3D::projectBottomToTop( const gp_Trsf &             bottom
     Handle(Geom_Surface) surface = BRep_Tool::Surface( topFace, loc );
     bool isPlanar = GeomLib_IsPlanarSurface( surface ).IsPlanar();
 
-    bool isFixed = false;
     set<const SMDS_MeshNode*> fixedNodes;
-    for ( int iAttemp = 0; !isFixed && iAttemp < 10; ++iAttemp )
-    {
-      TIDSortedElemSet faces;
-      for ( faceIt = topSMDS->GetElements(); faceIt->more(); )
-        faces.insert( faces.end(), faceIt->next() );
+    TIDSortedElemSet faces;
+    for ( faceIt = topSMDS->GetElements(); faceIt->more(); )
+      faces.insert( faces.end(), faceIt->next() );
 
+    bool isOk = false;
+    for ( int isCentroidal = 0; isCentroidal < 2; ++isCentroidal )
+    {
       SMESH_MeshEditor::SmoothMethod algo =
-        iAttemp ? SMESH_MeshEditor::CENTROIDAL : SMESH_MeshEditor::LAPLACIAN;
+        isCentroidal ? SMESH_MeshEditor::CENTROIDAL : SMESH_MeshEditor::LAPLACIAN;
+
+      int nbAttempts = isCentroidal ? 1 : 10;
+      for ( int iAttemp = 0; iAttemp < nbAttempts; ++iAttemp )
+      {
+        TIDSortedElemSet workFaces = faces;
 
-      // smoothing
-      editor.Smooth( faces, fixedNodes, algo, /*nbIterations=*/ 10,
-                     /*theTgtAspectRatio=*/1.0, /*the2D=*/!isPlanar);
+        // smoothing
+        editor.Smooth( workFaces, fixedNodes, algo, /*nbIterations=*/ 10,
+                       /*theTgtAspectRatio=*/1.0, /*the2D=*/!isPlanar);
 
-      isFixed = !topHelper.IsDistorted2D( topSM, /*checkUV=*/true );
+        if (( isOk = !topHelper.IsDistorted2D( topSM, /*checkUV=*/true )) &&
+            ( !isCentroidal ))
+          break;
+      }
     }
-    if ( !isFixed )
+    if ( !isOk )
       return toSM( error( TCom("Projection from face #") << botSM->GetId()
                           << " to face #" << topSM->GetId()
                           << " failed: inverted elements created"));
   }
 
+  TProjction2dAlgo::instance( this )->SetEventListener( topSM );
+
   return true;
 }
 
@@ -2356,6 +2423,9 @@ double StdMeshers_Prism_3D::getSweepTolerance( const Prism_3D::TPrismTopo& thePr
 
 bool StdMeshers_Prism_3D::isSimpleBottom( const Prism_3D::TPrismTopo& thePrism )
 {
+  if ( thePrism.myNbEdgesInWires.front() != 4 )
+    return false;
+
   // analyse angles between edges
   double nbConcaveAng = 0, nbConvexAng = 0;
   TopoDS_Face reverseBottom = TopoDS::Face( thePrism.myBottom.Reversed() ); // see initPrism()
@@ -2385,6 +2455,43 @@ bool StdMeshers_Prism_3D::isSimpleBottom( const Prism_3D::TPrismTopo& thePrism )
   return true;
 }
 
+//=======================================================================
+//function : allVerticalEdgesStraight
+//purpose  : Defines if all "vertical" EDGEs are straight
+//=======================================================================
+
+bool StdMeshers_Prism_3D::allVerticalEdgesStraight( const Prism_3D::TPrismTopo& thePrism )
+{
+  for ( size_t i = 0; i < thePrism.myWallQuads.size(); ++i )
+  {
+    const Prism_3D::TQuadList& quads = thePrism.myWallQuads[i];
+    Prism_3D::TQuadList::const_iterator quadIt = quads.begin();
+    TopoDS_Edge prevQuadEdge;
+    for ( ; quadIt != quads.end(); ++quadIt )
+    {
+      StdMeshers_FaceSidePtr rightSide = (*quadIt)->side[ QUAD_RIGHT_SIDE ];
+
+      if ( !prevQuadEdge.IsNull() &&
+           !SMESH_Algo::IsContinuous( rightSide->Edge( 0 ), prevQuadEdge ))
+        return false;
+
+      for ( int iE = 0; iE < rightSide->NbEdges(); ++iE )
+      {
+        const TopoDS_Edge & rightE = rightSide->Edge( iE );
+        if ( !SMESH_Algo::IsStraight( rightE, /*degenResult=*/true ))
+          return false;
+
+        if ( iE > 0 &&
+             !SMESH_Algo::IsContinuous( rightSide->Edge( iE-1 ), rightE ))
+          return false;
+
+        prevQuadEdge = rightE;
+      }
+    }
+  }
+  return true;
+}
+
 //=======================================================================
 //function : project2dMesh
 //purpose  : Project mesh faces from a source FACE of one prism (theSrcFace)
@@ -2411,6 +2518,8 @@ bool StdMeshers_Prism_3D::project2dMesh(const TopoDS_Face& theSrcFace,
   tgtSM->ComputeStateEngine       ( SMESH_subMesh::CHECK_COMPUTE_STATE );
   tgtSM->ComputeSubMeshStateEngine( SMESH_subMesh::CHECK_COMPUTE_STATE );
 
+  projector2D->SetEventListener( tgtSM );
+
   return ok;
 }
 
@@ -2490,26 +2599,34 @@ namespace // utils used by StdMeshers_Prism_3D::IsApplicable()
   struct EdgeWithNeighbors
   {
     TopoDS_Edge _edge;
-    int         _iL, _iR;
-    EdgeWithNeighbors(const TopoDS_Edge& E, int iE, int nbE, int shift = 0 ):
-      _edge( E ),
-      _iL( SMESH_MesherHelper::WrapIndex( iE-1, nbE ) + shift ),
-      _iR( SMESH_MesherHelper::WrapIndex( iE+1, nbE ) + shift )
+    int         _iBase;   /* index in a WIRE with non-base EDGEs excluded */
+    int         _iL, _iR; /* used to connect edges in a base FACE */
+    bool        _isBase;  /* is used in a base FACE */
+    EdgeWithNeighbors(const TopoDS_Edge& E, int iE, int nbE, int shift, bool isBase ):
+      _edge( E ), _iBase( iE + shift ),
+      _iL( SMESH_MesherHelper::WrapIndex( iE-1, Max( 1, nbE )) + shift ),
+      _iR( SMESH_MesherHelper::WrapIndex( iE+1, Max( 1, nbE )) + shift ),
+      _isBase( isBase )
     {
     }
     EdgeWithNeighbors() {}
+    bool IsInternal() const { return !_edge.IsNull() && _edge.Orientation() == TopAbs_INTERNAL; }
   };
-  struct PrismSide
+  // PrismSide contains all FACEs linking a bottom EDGE with a top one. 
+  struct PrismSide 
   {
-    TopoDS_Face                 _face;
-    TopTools_IndexedMapOfShape *_faces; // pointer because its copy constructor is private
-    TopoDS_Edge                 _topEdge;
-    vector< EdgeWithNeighbors >*_edges;
-    int                         _iBotEdge;
-    vector< bool >              _isCheckedEdge;
+    TopoDS_Face                 _face;    // a currently treated upper FACE
+    TopTools_IndexedMapOfShape *_faces;   // all FACEs (pointer because of a private copy constructor)
+    TopoDS_Edge                 _topEdge; // a current top EDGE
+    vector< EdgeWithNeighbors >*_edges;   // all EDGEs of _face
+    int                         _iBotEdge;       // index of _topEdge within _edges
+    vector< bool >              _isCheckedEdge;  // mark EDGEs whose two owner FACEs found
     int                         _nbCheckedEdges; // nb of EDGEs whose location is defined
-    PrismSide                  *_leftSide;
+    PrismSide                  *_leftSide;       // neighbor sides
     PrismSide                  *_rightSide;
+    bool                        _isInternal; // whether this side raises from an INTERNAL EDGE
+    void SetExcluded() { _leftSide = _rightSide = NULL; }
+    bool IsExcluded() const { return !_leftSide; }
     const TopoDS_Edge& Edge( int i ) const
     {
       return (*_edges)[ i ]._edge;
@@ -2520,67 +2637,138 @@ namespace // utils used by StdMeshers_Prism_3D::IsApplicable()
         if ( E.IsSame( Edge( i ))) return i;
       return -1;
     }
-    bool IsSideFace( const TopoDS_Shape& face ) const
+    bool IsSideFace( const TopoDS_Shape& face, const bool checkNeighbors ) const
     {
       if ( _faces->Contains( face )) // avoid returning true for a prism top FACE
         return ( !_face.IsNull() || !( face.IsSame( _faces->FindKey( _faces->Extent() ))));
+
+      if ( checkNeighbors )
+        return (( _leftSide  && _leftSide->IsSideFace ( face, false )) ||
+                ( _rightSide && _rightSide->IsSideFace( face, false )));
+
       return false;
     }
   };
+  //--------------------------------------------------------------------------------
+  /*!
+   * \brief Return another faces sharing an edge
+   */
+  const TopoDS_Face & getAnotherFace( const TopoDS_Face& face,
+                                      const TopoDS_Edge& edge,
+                                      TopTools_IndexedDataMapOfShapeListOfShape& facesOfEdge)
+  {
+    TopTools_ListIteratorOfListOfShape faceIt( facesOfEdge.FindFromKey( edge ));
+    for ( ; faceIt.More(); faceIt.Next() )
+      if ( !face.IsSame( faceIt.Value() ))
+        return TopoDS::Face( faceIt.Value() );
+    return face;
+  }
+
   //--------------------------------------------------------------------------------
   /*!
    * \brief Return ordered edges of a face
    */
-  bool getEdges( const TopoDS_Face&            face,
-                 vector< EdgeWithNeighbors > & edges,
-                 const bool                    noHolesAllowed)
+  bool getEdges( const TopoDS_Face&                         face,
+                 vector< EdgeWithNeighbors > &              edges,
+                 TopTools_IndexedDataMapOfShapeListOfShape& facesOfEdge,
+                 const bool                                 noHolesAllowed)
   {
+    TopoDS_Face f = face;
+    if ( f.Orientation() != TopAbs_FORWARD &&
+         f.Orientation() != TopAbs_REVERSED )
+      f.Orientation( TopAbs_FORWARD );
     list< TopoDS_Edge > ee;
     list< int >         nbEdgesInWires;
-    int nbW = SMESH_Block::GetOrderedEdges( face, ee, nbEdgesInWires );
+    int nbW = SMESH_Block::GetOrderedEdges( f, ee, nbEdgesInWires );
     if ( nbW > 1 && noHolesAllowed )
       return false;
 
-    int iE, nbTot = 0;
-    list< TopoDS_Edge >::iterator e = ee.begin();
-    list< int >::iterator       nbE = nbEdgesInWires.begin();
+    int iE, nbTot = 0, nbBase, iBase;
+    list< TopoDS_Edge >::iterator   e = ee.begin();
+    list< int         >::iterator nbE = nbEdgesInWires.begin();
     for ( ; nbE != nbEdgesInWires.end(); ++nbE )
       for ( iE = 0; iE < *nbE; ++e, ++iE )
-        if ( SMESH_Algo::isDegenerated( *e ))
+        if ( SMESH_Algo::isDegenerated( *e )) // degenerated EDGE is never used
         {
           e = --ee.erase( e );
           --(*nbE);
           --iE;
         }
-        else
-        {
-          e->Orientation( TopAbs_FORWARD ); // for operator==() to work
-        }
 
+    vector<int> isBase;
     edges.clear();
     e = ee.begin();
     for ( nbE = nbEdgesInWires.begin(); nbE != nbEdgesInWires.end(); ++nbE )
     {
-      for ( iE = 0; iE < *nbE; ++e, ++iE )
-        edges.push_back( EdgeWithNeighbors( *e, iE, *nbE, nbTot ));
-      nbTot += *nbE;
+      nbBase = 0;
+      isBase.resize( *nbE );
+      list< TopoDS_Edge >::iterator eIt = e;
+      for ( iE = 0; iE < *nbE; ++eIt, ++iE )
+      {
+        isBase[ iE ] = ( getAnotherFace( face, *eIt, facesOfEdge ) != face );
+        nbBase += isBase[ iE ];
+      }
+      for ( iBase = 0, iE = 0; iE < *nbE; ++e, ++iE )
+      {
+        edges.push_back( EdgeWithNeighbors( *e, iBase, nbBase, nbTot, isBase[ iE ] ));
+        iBase += isBase[ iE ];
+      }
+      nbTot += nbBase;
+    }
+    if ( nbTot == 0 )
+      return false;
+
+    // IPAL53099. Set correct neighbors to INTERNAL EDGEs, which can be connected to
+    // EDGEs of the outer WIRE but this fact can't be detected by their order.
+    if ( nbW > 1 )
+    {
+      int iFirst = 0, iLast;
+      for ( nbE = nbEdgesInWires.begin(); nbE != nbEdgesInWires.end(); ++nbE )
+      {
+        iLast = iFirst + *nbE - 1;
+        TopoDS_Vertex vv[2] = { SMESH_MesherHelper::IthVertex( 0, edges[ iFirst ]._edge ),
+                                SMESH_MesherHelper::IthVertex( 1, edges[ iLast  ]._edge ) };
+        bool isConnectOk = ( vv[0].IsSame( vv[1] ));
+        if ( !isConnectOk )
+        {
+          edges[ iFirst ]._iL = edges[ iFirst ]._iBase; // connect to self
+          edges[ iLast  ]._iR = edges[ iLast ]._iBase;
+
+          // look for an EDGE of the outer WIREs connected to vv
+          TopoDS_Vertex v0, v1;
+          for ( iE = 0; iE < iFirst; ++iE )
+          {
+            v0 = SMESH_MesherHelper::IthVertex( 0, edges[ iE ]._edge );
+            v1 = SMESH_MesherHelper::IthVertex( 1, edges[ iE ]._edge );
+            if ( vv[0].IsSame( v0 ) || vv[0].IsSame( v1 ))
+              edges[ iFirst ]._iL = edges[ iE ]._iBase;
+            if ( vv[1].IsSame( v0 ) || vv[1].IsSame( v1 ))
+              edges[ iLast  ]._iR = edges[ iE ]._iBase;
+          }
+        }
+        iFirst += *nbE;
+      }
     }
     return edges.size();
   }
+  
   //--------------------------------------------------------------------------------
   /*!
-   * \brief Return another faces sharing an edge
+   * \brief Return number of faces sharing given edges
    */
-  const TopoDS_Shape & getAnotherFace( const TopoDS_Face& face,
-                                       const TopoDS_Edge& edge,
-                                       TopTools_IndexedDataMapOfShapeListOfShape& facesOfEdge)
-  {
-    TopTools_ListIteratorOfListOfShape faceIt( facesOfEdge.FindFromKey( edge ));
-    for ( ; faceIt.More(); faceIt.Next() )
-      if ( !face.IsSame( faceIt.Value() ))
-        return faceIt.Value();
-    return face;
-  }
+  // int nbAdjacentFaces( const std::vector< EdgeWithNeighbors >&          edges,
+  //                      const TopTools_IndexedDataMapOfShapeListOfShape& facesOfEdge )
+  // {
+  //   TopTools_MapOfShape adjFaces;
+
+  //   for ( size_t i = 0; i < edges.size(); ++i )
+  //   {
+  //     TopTools_ListIteratorOfListOfShape faceIt( facesOfEdge.FindFromKey( edges[i]._edge ));
+  //     for ( ; faceIt.More(); faceIt.Next() )
+  //       adjFaces.Add( faceIt.Value() );
+  //   }
+  //   return adjFaces.Extent();
+  // }
 }
 
 //================================================================================
@@ -2603,10 +2791,10 @@ bool StdMeshers_Prism_3D::IsApplicable(const TopoDS_Shape & shape, bool toCheckA
     // check nb shells
     TopoDS_Shape shell;
     TopExp_Explorer shExp( sExp.Current(), TopAbs_SHELL );
-    if ( shExp.More() ) {
+    while ( shExp.More() ) {
       shell = shExp.Current();
       shExp.Next();
-      if ( shExp.More() )
+      if ( shExp.More() && BRep_Tool::IsClosed( shExp.Current() ))
         shell.Nullify();
     }
     if ( shell.IsNull() ) {
@@ -2615,7 +2803,7 @@ bool StdMeshers_Prism_3D::IsApplicable(const TopoDS_Shape & shape, bool toCheckA
     }
     // get all faces
     TopTools_IndexedMapOfShape allFaces;
-    TopExp::MapShapes( shell, TopAbs_FACE, allFaces );
+    TopExp::MapShapes( sExp.Current(), TopAbs_FACE, allFaces );
     if ( allFaces.Extent() < 3 ) {
       if ( toCheckAll ) return false;
       continue;
@@ -2632,7 +2820,7 @@ bool StdMeshers_Prism_3D::IsApplicable(const TopoDS_Shape & shape, bool toCheckA
       }
     }
 #ifdef _DEBUG_
-    TopTools_IndexedMapOfShape allShapes;
+    TopTools_IndexedMapOfShape allShapes; // usage: allShapes.FindIndex( s )
     TopExp::MapShapes( shape, allShapes );
 #endif
 
@@ -2646,32 +2834,45 @@ bool StdMeshers_Prism_3D::IsApplicable(const TopoDS_Shape & shape, bool toCheckA
 
     typedef vector< EdgeWithNeighbors > TEdgeWithNeighborsVec;
     vector< TEdgeWithNeighborsVec > faceEdgesVec( allFaces.Extent() + 1 );
-    TopTools_IndexedMapOfShape* facesOfSide = new TopTools_IndexedMapOfShape[ faceEdgesVec.size() ];
+    const size_t nbEdgesMax = facesOfEdge.Extent() * 2; // there can be seam EDGEs
+    TopTools_IndexedMapOfShape* facesOfSide = new TopTools_IndexedMapOfShape[ nbEdgesMax ];
     SMESHUtils::ArrayDeleter<TopTools_IndexedMapOfShape> delFacesOfSide( facesOfSide );
 
     // try to use each face as a bottom one
     bool prismDetected = false;
+    vector< PrismSide > sides;
     for ( int iF = 1; iF < allFaces.Extent() && !prismDetected; ++iF )
     {
       const TopoDS_Face& botF = TopoDS::Face( allFaces( iF ));
 
       TEdgeWithNeighborsVec& botEdges = faceEdgesVec[ iF ];
       if ( botEdges.empty() )
-        if ( !getEdges( botF, botEdges, /*noHoles=*/false ))
+        if ( !getEdges( botF, botEdges, facesOfEdge, /*noHoles=*/false ))
           break;
-      if ( allFaces.Extent()-1 <= (int) botEdges.size() )
+
+      int nbBase = 0;
+      for ( size_t iS = 0; iS < botEdges.size(); ++iS )
+        nbBase += botEdges[ iS ]._isBase;
+
+      if ( allFaces.Extent()-1 <= nbBase )
         continue; // all faces are adjacent to botF - no top FACE
 
       // init data of side FACEs
-      vector< PrismSide > sides( botEdges.size() );
-      for ( int iS = 0; iS < botEdges.size(); ++iS )
+      sides.clear();
+      sides.resize( nbBase );
+      size_t iS = 0;
+      for ( size_t iE = 0; iE < botEdges.size(); ++iE )
       {
-        sides[ iS ]._topEdge = botEdges[ iS ]._edge;
-        sides[ iS ]._face    = botF;
-        sides[ iS ]._leftSide  = & sides[ botEdges[ iS ]._iR ];
-        sides[ iS ]._rightSide = & sides[ botEdges[ iS ]._iL ];
-        sides[ iS ]._faces = & facesOfSide[ iS ];
+        if ( !botEdges[ iE ]._isBase )
+          continue;
+        sides[ iS ]._topEdge    = botEdges[ iE ]._edge;
+        sides[ iS ]._face       = botF;
+        sides[ iS ]._leftSide   = & sides[ botEdges[ iE ]._iR ];
+        sides[ iS ]._rightSide  = & sides[ botEdges[ iE ]._iL ];
+        sides[ iS ]._isInternal = botEdges[ iE ].IsInternal();
+        sides[ iS ]._faces      = & facesOfSide[ iS ];
         sides[ iS ]._faces->Clear();
+        ++iS;
       }
 
       bool isOK = true; // ok for a current botF
@@ -2699,8 +2900,9 @@ bool StdMeshers_Prism_3D::IsApplicable(const TopoDS_Shape & shape, bool toCheckA
                 if ( side._isCheckedEdge[ iE ] ) continue;
                 const TopoDS_Edge&      vertE = side.Edge( iE );
                 const TopoDS_Shape& neighborF = getAnotherFace( side._face, vertE, facesOfEdge );
-                bool isEdgeShared = adjSide->IsSideFace( neighborF );
-                if ( isEdgeShared )
+                bool isEdgeShared = (( adjSide->IsSideFace( neighborF, side._isInternal )) ||
+                                     ( adjSide == &side && neighborF.IsSame( side._face )) );
+                if ( isEdgeShared ) // vertE is shared with adjSide
                 {
                   isAdvanced = true;
                   side._isCheckedEdge[ iE ] = true;
@@ -2741,20 +2943,19 @@ bool StdMeshers_Prism_3D::IsApplicable(const TopoDS_Shape & shape, bool toCheckA
             {
               stop = true;
             }
-            else if ( side._leftSide != & side ) // not closed side face
+            else if ( side._leftSide != & side && // not closed side face
+                      side._leftSide->_faces->Contains( f ))
             {
-              if ( side._leftSide->_faces->Contains( f ))
-              {
-                stop = true; // probably f is the prism top face
-                side._leftSide->_face.Nullify();
-                side._leftSide->_topEdge.Nullify();
-              }
-              if ( side._rightSide->_faces->Contains( f ))
-              {
-                stop = true; // probably f is the prism top face
-                side._rightSide->_face.Nullify();
-                side._rightSide->_topEdge.Nullify();
-              }
+              stop = true; // probably f is the prism top face
+              side._leftSide->_face.Nullify();
+              side._leftSide->_topEdge.Nullify();
+            }
+            else if ( side._rightSide != & side &&
+                      side._rightSide->_faces->Contains( f ))
+            {
+              stop = true; // probably f is the prism top face
+              side._rightSide->_face.Nullify();
+              side._rightSide->_topEdge.Nullify();
             }
             if ( stop )
             {
@@ -2766,7 +2967,7 @@ bool StdMeshers_Prism_3D::IsApplicable(const TopoDS_Shape & shape, bool toCheckA
             int faceID  = allFaces.FindIndex( side._face );
             side._edges = & faceEdgesVec[ faceID ];
             if ( side._edges->empty() )
-              if ( !getEdges( side._face, * side._edges, /*noHoles=*/true ))
+              if ( !getEdges( side._face, * side._edges, facesOfEdge, /*noHoles=*/true ))
                 break;
             const int nbE = side._edges->size();
             if ( nbE >= 4 )
@@ -2779,6 +2980,10 @@ bool StdMeshers_Prism_3D::IsApplicable(const TopoDS_Shape & shape, bool toCheckA
               side._isCheckedEdge[ side._iBotEdge ] = true;
               side._nbCheckedEdges = 1; // bottom EDGE is known
             }
+            else // probably a triangular top face found
+            {
+              side._face.Nullify();
+            }
             side._topEdge.Nullify();
             isOK = ( !side._edges->empty() || side._faces->Extent() > 1 );
 
@@ -2811,7 +3016,7 @@ bool StdMeshers_Prism_3D::IsApplicable(const TopoDS_Shape & shape, bool toCheckA
           const TopoDS_Shape& topFace = sides[0]._faces->FindKey( nbFaces );
           size_t iS;
           for ( iS = 1; iS < sides.size(); ++iS )
-            if ( !sides[ iS ]._faces->Contains( topFace ))
+            if ( ! sides[ iS ]._faces->Contains( topFace ))
               break;
           prismDetected = ( iS == sides.size() );
         }
@@ -2937,7 +3142,6 @@ bool StdMeshers_Prism_3D::initPrism(Prism_3D::TPrismTopo& thePrism,
   list< SMESH_subMesh* > meshedSubMesh;
   int nbFaces = 0;
   //
-  SMESH_subMesh* anyFaceSM = 0;
   SMESH_subMeshIteratorPtr smIt = mainSubMesh->getDependsOnIterator(false,true);
   while ( smIt->more() )
   {
@@ -2946,7 +3150,6 @@ bool StdMeshers_Prism_3D::initPrism(Prism_3D::TPrismTopo& thePrism,
     if      ( face.ShapeType() > TopAbs_FACE ) break;
     else if ( face.ShapeType() < TopAbs_FACE ) continue;
     nbFaces++;
-    anyFaceSM = sm;
 
     // is quadrangle FACE?
     list< TopoDS_Edge > orderedEdges;
@@ -3076,7 +3279,7 @@ bool StdMeshers_Prism_3D::initPrism(Prism_3D::TPrismTopo& thePrism,
 
   // find vertex 000 - the one with smallest coordinates (for easy DEBUG :-)
   TopoDS_Vertex V000;
-  double minVal = DBL_MAX, minX, val;
+  double minVal = DBL_MAX, minX = 0, val;
   for ( TopExp_Explorer exp( botSM->GetSubShape(), TopAbs_VERTEX );
         exp.More(); exp.Next() )
   {
@@ -3115,9 +3318,11 @@ bool StdMeshers_Prism_3D::initPrism(Prism_3D::TPrismTopo& thePrism,
                     "Non-quadrilateral faces are not opposite"));
 
     // check that the found top and bottom FACEs are opposite
+    TopTools_IndexedMapOfShape topEdgesMap( thePrism.myBottomEdges.size() );
+    TopExp::MapShapes( thePrism.myTop, topEdgesMap );
     list< TopoDS_Edge >::iterator edge = thePrism.myBottomEdges.begin();
     for ( ; edge != thePrism.myBottomEdges.end(); ++edge )
-      if ( myHelper->IsSubShape( *edge, thePrism.myTop ))
+      if ( topEdgesMap.Contains( *edge ))
         return toSM( error
                      (notQuadGeomSubMesh.empty() ? COMPERR_BAD_INPUT_MESH : COMPERR_BAD_SHAPE,
                       "Non-quadrilateral faces are not opposite"));
@@ -3219,6 +3424,9 @@ bool StdMeshers_PrismAsBlock::Init(SMESH_MesherHelper*         helper,
       if ( !myHelper->LoadNodeColumns( faceColumns, (*quad)->face, quadBot, meshDS ))
         return error(COMPERR_BAD_INPUT_MESH, TCom("Can't find regular quadrangle mesh ")
                      << "on a side face #" << MeshDS()->ShapeToIndex( (*quad)->face ));
+
+      if ( !faceColumns.empty() && (int)faceColumns.begin()->second.size() != VerticalSize() )
+        return error(COMPERR_BAD_INPUT_MESH, "Different 'vertical' discretization");
     }
     // edge columns
     int id = MeshDS()->ShapeToIndex( *edgeIt );
@@ -3285,7 +3493,7 @@ bool StdMeshers_PrismAsBlock::Init(SMESH_MesherHelper*         helper,
       if ( nbUnitePerEdge[ iE ] < 0 )
         continue;
       // look for already united faces
-      for ( int i = iE; i < iE + nbExraFaces; ++i )
+      for ( size_t i = iE; i < iE + nbExraFaces; ++i )
       {
         if ( nbUnitePerEdge[ i ] > 0 ) // a side including nbUnitePerEdge[i]+1 edge
           nbExraFaces += nbUnitePerEdge[ i ];
@@ -3324,7 +3532,7 @@ bool StdMeshers_PrismAsBlock::Init(SMESH_MesherHelper*         helper,
     else if ( nbExraFaces > 1 ) // unite
     {
       double u0 = 0, sumLen = 0;
-      for ( int i = iE; i < iE + nbExraFaces; ++i )
+      for ( size_t i = iE; i < iE + nbExraFaces; ++i )
         sumLen += edgeLength[ i ];
 
       vector< TSideFace* >        components( nbExraFaces );
@@ -3568,7 +3776,7 @@ bool StdMeshers_PrismAsBlock::GetLayersTransformation(vector<gp_Trsf> &
   double tol2;
   {
     Bnd_B3d bndBox;
-    for ( int i = 0; i < columns.size(); ++i )
+    for ( size_t i = 0; i < columns.size(); ++i )
       bndBox.Add( gpXYZ( columns[i]->front() ));
     tol2 = bndBox.SquareExtent() * 1e-5;
   }
@@ -3591,7 +3799,7 @@ bool StdMeshers_PrismAsBlock::GetLayersTransformation(vector<gp_Trsf> &
     //t.SetScaleFactor( distZ/dist0 ); - it does not work properly, wrong base point
 
     // check a transformation
-    for ( int i = 0; i < columns.size(); ++i )
+    for ( size_t i = 0; i < columns.size(); ++i )
     {
       gp_Pnt p0 = gpXYZ( (*columns[i])[0] );
       gp_Pnt pz = gpXYZ( (*columns[i])[z] );
@@ -3795,7 +4003,7 @@ StdMeshers_PrismAsBlock::TSideFace::TSideFace( const TSideFace& other ):
   myComponents       ( other.myComponents.size() ),
   myHelper           ( *other.myHelper.GetMesh() )
 {
-  for (int i = 0 ; i < myComponents.size(); ++i )
+  for ( size_t i = 0 ; i < myComponents.size(); ++i )
     myComponents[ i ] = new TSideFace( *other.myComponents[ i ]);
 }
 
@@ -3807,7 +4015,7 @@ StdMeshers_PrismAsBlock::TSideFace::TSideFace( const TSideFace& other ):
 
 StdMeshers_PrismAsBlock::TSideFace::~TSideFace()
 {
-  for (int i = 0 ; i < myComponents.size(); ++i )
+  for ( size_t i = 0 ; i < myComponents.size(); ++i )
     if ( myComponents[ i ] )
       delete myComponents[ i ];
 }
@@ -3907,7 +4115,7 @@ StdMeshers_PrismAsBlock::TSideFace::GetComponent(const double U,double & localU)
   if ( myComponents.empty() )
     return const_cast<TSideFace*>( this );
 
-  int i;
+  size_t i;
   for ( i = 0; i < myComponents.size(); ++i )
     if ( U < myParams[ i ].second )
       break;
@@ -4346,9 +4554,9 @@ gp_Pnt StdMeshers_PrismAsBlock::TVerticalEdgeAdaptor::Value(const Standard_Real
 void StdMeshers_PrismAsBlock::TVerticalEdgeAdaptor::dumpNodes(int nbNodes) const
 {
 #ifdef _DEBUG_
-  for ( int i = 0; i < nbNodes && i < myNodeColumn->size(); ++i )
+  for ( int i = 0; i < nbNodes && i < (int)myNodeColumn->size(); ++i )
     cout << (*myNodeColumn)[i]->GetID() << " ";
-  if ( nbNodes < myNodeColumn->size() )
+  if ( nbNodes < (int) myNodeColumn->size() )
     cout << myNodeColumn->back()->GetID();
 #endif
 }
@@ -4584,7 +4792,7 @@ void StdMeshers_Sweeper::applyBoundaryError(const vector< gp_XYZ >& bndPoints,
 
 //================================================================================
 /*!
- * \brief Creates internal nodes of the prism
+ * \brief Create internal nodes of the prism
  */
 //================================================================================
 
@@ -4746,8 +4954,8 @@ bool StdMeshers_Sweeper::ComputeNodes( SMESH_MesherHelper& helper,
     }
   }
 
-  //centerIntErrorIsSmall = true;
-  //bndErrorIsSmall = true;
+  centerIntErrorIsSmall = true; // 3D_mesh_Extrusion_00/A3
+  bndErrorIsSmall = true;
   if ( !centerIntErrorIsSmall )
   {
     // Compensate the central error; continue adding projection
@@ -4864,3 +5072,257 @@ bool StdMeshers_Sweeper::ComputeNodes( SMESH_MesherHelper& helper,
 
   return true;
 }
+
+//================================================================================
+/*!
+ * \brief Check if all nodes of each layers have same logical Z
+ */
+//================================================================================
+
+bool StdMeshers_Sweeper::CheckSameZ()
+{
+  myZColumns.resize( myBndColumns.size() );
+  fillZColumn( myZColumns[0], *myBndColumns[0] );
+
+  bool sameZ = true;
+  const double tol = 0.1 * 1./ myBndColumns[0]->size();
+
+  // check columns based on VERTEXes
+
+  vector< int > vertexIndex;
+  vertexIndex.push_back( 0 );
+  for ( size_t iC = 1; iC < myBndColumns.size() &&  sameZ; ++iC )
+  {
+    if ( myBndColumns[iC]->front()->GetPosition()->GetDim() > 0 )
+      continue; // not on VERTEX
+
+    vertexIndex.push_back( iC );
+    fillZColumn( myZColumns[iC], *myBndColumns[iC] );
+
+    for ( size_t iZ = 0; iZ < myZColumns[0].size() &&  sameZ; ++iZ )
+      sameZ = ( Abs( myZColumns[0][iZ] - myZColumns[iC][iZ]) < tol );
+  }
+
+  // check columns based on EDGEs, one per EDGE
+
+  for ( size_t i = 1; i < vertexIndex.size() &&  sameZ; ++i )
+  {
+    if ( vertexIndex[i] - vertexIndex[i-1] < 2 )
+      continue;
+
+    int iC = ( vertexIndex[i] + vertexIndex[i-1] ) / 2;
+    fillZColumn( myZColumns[iC], *myBndColumns[iC] );
+
+    for ( size_t iZ = 0; iZ < myZColumns[0].size() &&  sameZ; ++iZ )
+      sameZ = ( Abs( myZColumns[0][iZ] - myZColumns[iC][iZ]) < tol );
+  }
+
+  if ( sameZ )
+  {
+    myZColumns.resize(1);
+  }
+  else
+  {
+    for ( size_t iC = 1; iC < myBndColumns.size(); ++iC )
+      fillZColumn( myZColumns[iC], *myBndColumns[iC] );
+  }
+
+  return sameZ;
+}
+
+//================================================================================
+/*!
+ * \brief Create internal nodes of the prism all located on straight lines with
+ *        the same distribution along the lines.
+ */
+//================================================================================
+
+bool StdMeshers_Sweeper::ComputeNodesOnStraightSameZ( SMESH_MesherHelper& helper )
+{
+  TZColumn& z = myZColumns[0];
+
+  for ( size_t i = 0; i < myIntColumns.size(); ++i )
+  {
+    TNodeColumn& nodes = *myIntColumns[i];
+    SMESH_NodeXYZ n0( nodes[0] ), n1( nodes.back() );
+
+    for ( size_t iZ = 0; iZ < z.size(); ++iZ )
+    {
+      gp_XYZ p = n0 * ( 1 - z[iZ] ) + n1 * z[iZ];
+      nodes[ iZ+1 ] = helper.AddNode( p.X(), p.Y(), p.Z() );
+    }
+  }
+
+  return true;
+}
+
+//================================================================================
+/*!
+ * \brief Create internal nodes of the prism all located on straight lines with
+ *        different distributions along the lines.
+ */
+//================================================================================
+
+bool StdMeshers_Sweeper::ComputeNodesOnStraight( SMESH_MesherHelper& helper,
+                                                 const TopoDS_Face&  botFace,
+                                                 const TopoDS_Face&  topFace )
+{
+  // get data to create a Morph
+  UVPtStructVec botUV( myBndColumns.size() + 1 );
+  UVPtStructVec topUV( myBndColumns.size() + 1 );
+  for ( size_t i = 0; i < myBndColumns.size(); ++i )
+  {
+    TNodeColumn& nodes = *myBndColumns[i];
+    botUV[i].node = nodes[0];
+    botUV[i].SetUV( helper.GetNodeUV( botFace, nodes[0] ));
+    topUV[i].node = nodes.back();
+    topUV[i].SetUV( helper.GetNodeUV( topFace, nodes.back() ));
+    botUV[i].node->setIsMarked( true );
+  }
+  botUV.back() = botUV[0];
+  topUV.back() = topUV[0];
+
+  TopoDS_Edge dummyE;
+  TSideVector botWires( 1, StdMeshers_FaceSide::New( botUV, botFace, dummyE, helper.GetMesh() ));
+  TSideVector topWires( 1, StdMeshers_FaceSide::New( topUV, topFace, dummyE, helper.GetMesh() ));
+
+  // use Morph to make delauney mesh on the FACEs. Locating of a node within a
+  // delauney triangle will be used to get a weighted Z.
+  NSProjUtils::Morph botDelauney( botWires );
+  NSProjUtils::Morph topDelauney( topWires );
+
+  if ( helper.GetIsQuadratic() )
+  {
+    // mark all medium nodes of faces on botFace to avoid their treating
+    SMESHDS_SubMesh* smDS = helper.GetMeshDS()->MeshElements( botFace );
+    SMDS_ElemIteratorPtr eIt = smDS->GetElements();
+    while ( eIt->more() )
+    {
+      const SMDS_MeshElement* e = eIt->next();
+      for ( int i = e->NbCornerNodes(), nb = e->NbNodes(); i < nb; ++i )
+        e->GetNode( i )->setIsMarked( true );
+    }
+  }
+
+  // map to get a node column by a bottom node
+  TColStd_DataMapOfIntegerInteger iNode2iCol( myIntColumns.size() );
+
+  // un-mark nodes to treat (internal bottom nodes); later we will mark treated nodes
+  for ( size_t i = 0; i < myIntColumns.size(); ++i )
+  {
+    const SMDS_MeshNode* botNode = myIntColumns[i]->front();
+    botNode->setIsMarked( false );
+    iNode2iCol.Bind( botNode->GetID(), i );
+  }
+
+  const int botFaceID = helper.GetMesh()->GetSubMesh( botFace )->GetId();
+  const SMDS_MeshNode     *botNode, *topNode;
+  const BRepMesh_Triangle *botTria, *topTria;
+  double botBC[3], topBC[3]; // barycentric coordinates
+  int    botTriaNodes[3], topTriaNodes[3];
+  bool   checkUV = true;
+
+  // a queue of bottom nodes with starting delauney triangles
+  NSProjUtils::Morph::TNodeTriaList botNoTriQueue;
+
+  size_t iBndN = 1; // index of a bottom boundary node
+  int nbNodesToProcess = myIntColumns.size();
+  while ( nbNodesToProcess > 0 )
+  {
+    while ( !botNoTriQueue.empty() ) // treat all nodes in the queue
+    {
+      botNode = botNoTriQueue.front().first;
+      botTria = botNoTriQueue.front().second;
+      botNoTriQueue.pop_front();
+      if ( botNode->isMarked() )
+        continue;
+      --nbNodesToProcess;
+      botNode->setIsMarked( true );
+
+      TNodeColumn* column = myIntColumns[ iNode2iCol( botNode->GetID() )];
+
+      // find a delauney triangle containing the botNode
+      gp_XY botUV = helper.GetNodeUV( botFace, botNode, NULL, &checkUV );
+      botUV  *= botDelauney.GetScale();
+      botTria = botDelauney.FindTriangle( botUV, botTria, botBC, botTriaNodes );
+      if ( !botTria )
+        return false;
+
+      // find a delauney triangle containing the topNode
+      topNode = column->back();
+      gp_XY topUV = helper.GetNodeUV( topFace, topNode, NULL, &checkUV );
+      topUV *= topDelauney.GetScale();
+      // get a starting triangle basing on that top and bot boundary nodes have same index
+      topTria = topDelauney.GetTriangleNear( botTriaNodes[0] );
+      topTria = topDelauney.FindTriangle( topUV, topTria, topBC, topTriaNodes );
+      if ( !topTria )
+        return false;
+
+      // create nodes along a line
+      SMESH_NodeXYZ botP( botNode ), topP( topNode);
+      for ( size_t iZ = 0; iZ < myZColumns[0].size(); ++iZ )
+      {
+        // use barycentric coordinates as weight of Z of boundary columns
+        double botZ = 0, topZ = 0;
+        for ( int i = 0; i < 3; ++i )
+        {
+          botZ += botBC[i] * myZColumns[ botTriaNodes[i]-1 ][ iZ ];
+          topZ += topBC[i] * myZColumns[ topTriaNodes[i]-1 ][ iZ ];
+        }
+        double rZ = double( iZ + 1 ) / ( myZColumns[0].size() + 1 );
+        double z = botZ * ( 1 - rZ ) + topZ * rZ;
+        gp_XYZ p = botP * ( 1 - z  ) + topP * z;
+        (*column)[ iZ+1 ] = helper.AddNode( p.X(), p.Y(), p.Z() );
+      }
+
+      // add neighbor nodes to the queue
+      botDelauney.AddCloseNodes( botNode, botTria, botFaceID, botNoTriQueue );
+    }
+
+    if ( nbNodesToProcess > 0 ) // fill the queue
+    {
+      // assure that all bot nodes are visited
+      for ( ; iBndN-1 < myBndColumns.size() &&  botNoTriQueue.empty();  ++iBndN )
+      {
+        botTria = botDelauney.GetTriangleNear( iBndN );
+        const SMDS_MeshNode*  bndNode = botDelauney.GetBndNodes()[ iBndN ];
+        botDelauney.AddCloseNodes( bndNode, botTria, botFaceID, botNoTriQueue );
+      }
+      if ( botNoTriQueue.empty() )
+      {
+        for ( size_t i = 0; i < myIntColumns.size(); ++i )
+        {
+          botNode = myIntColumns[i]->front();
+          if ( !botNode->isMarked() )
+            botNoTriQueue.push_back( make_pair( botNode, botTria ));
+        }
+      }
+    }
+  }
+
+  return true;
+}
+
+//================================================================================
+/*!
+ * \brief Compute Z of nodes of a straight column
+ */
+//================================================================================
+
+void StdMeshers_Sweeper::fillZColumn( TZColumn&    zColumn,
+                                      TNodeColumn& nodes )
+{
+  if ( zColumn.size() == nodes.size() - 2 )
+    return;
+
+  gp_Pnt p0 = SMESH_NodeXYZ( nodes[0] );
+  gp_Vec line( p0, SMESH_NodeXYZ( nodes.back() ));
+  double len2 = line.SquareMagnitude();
+
+  zColumn.resize( nodes.size() - 2 );
+  for ( size_t i = 0; i < zColumn.size(); ++i )
+  {
+    gp_Vec vec( p0, SMESH_NodeXYZ( nodes[ i+1] ));
+    zColumn[i] = ( line * vec ) / len2; // param [0,1] on the line
+  }
+}