IPAL52499: Prismatic mesh is not computed on a prismatic shape

src/SMESHGUI/SMESHGUI.cxx

@@ -630,7 +630,7 @@ namespace
       };
       // is typeMsg complete? (compilation failure mains that enum SMDSAbs_EntityType changed)
       const int nbTypes = sizeof( typeMsg ) / sizeof( const char* );
-      int _assert[( nbTypes == SMESH::Entity_Last ) ? 2 : -1 ]; _assert[0]=_assert[1];
+      int _assert[( nbTypes == SMESH::Entity_Last ) ? 2 : -1 ]; _assert[0]=_assert[1]=0;
 
       QString andStr = " " + QObject::tr("SMESH_AND") + " ", comma(", ");
       for ( size_t iType = 0; iType < presentNotSupported.size(); ++iType ) {

src/SMESH_I/SMESH_Filter_i.cxx

@@ -4109,7 +4109,7 @@ static const char** getFunctNames()
 #ifdef _DEBUG_
   // check if functName is complete, compilation failure means that enum FunctorType changed
   const int nbFunctors = sizeof(functName) / sizeof(const char*);
-  int _assert[( nbFunctors == SMESH::FT_Undefined + 1 ) ? 2 : -1 ]; _assert[0]=_assert[1];
+  int _assert[( nbFunctors == SMESH::FT_Undefined + 1 ) ? 2 : -1 ]; _assert[0]=_assert[1]=0;
 #endif
 
   return functName;

src/SMESH_I/SMESH_Gen_i.cxx

@@ -2560,7 +2560,7 @@ SMESH_Gen_i::ConcatenateCommon(const SMESH::ListOfIDSources& theMeshesArray,
         const char* typeNames[] = { "All","Nodes","Edges","Faces","Volumes","0DElems","Balls" };
         { // check of typeNames: compilation failure mains that NB_ELEMENT_TYPES changed:
           const int nbNames = sizeof(typeNames) / sizeof(const char*);
-          int _assert[( nbNames == SMESH::NB_ELEMENT_TYPES ) ? 2 : -1 ]; _assert[0]=_assert[1];
+          int _assert[( nbNames == SMESH::NB_ELEMENT_TYPES ) ? 2 : -1 ]; _assert[0]=_assert[1]=0;
         }
         string groupName = "Gr";
         SALOMEDS::SObject_wrap aMeshSObj = ObjectToSObject( myCurrentStudy, theMeshesArray[i] );

src/SMESH_SWIG/smeshBuilder.py

@@ -2687,9 +2687,14 @@ class Mesh:
 
     ## Return an element based on all given nodes.
     #  @ingroup l1_meshinfo
-    def FindElementByNodes(self,nodes):
+    def FindElementByNodes(self, nodes):
         return self.mesh.FindElementByNodes(nodes)
 
+    ## Return elements including all given nodes.
+    #  @ingroup l1_meshinfo
+    def GetElementsByNodes(self, nodes, elemType=SMESH.ALL):
+        return self.mesh.GetElementsByNodes( nodes, elemType )
+
     ## Return true if the given element is a polygon
     #  @ingroup l1_meshinfo
     def IsPoly(self, id):

src/StdMeshers/StdMeshers_FaceSide.cxx

@@ -89,7 +89,7 @@ StdMeshers_FaceSide::StdMeshers_FaceSide(const TopoDS_Face&   theFace,
 //================================================================================
 
 StdMeshers_FaceSide::StdMeshers_FaceSide(const TopoDS_Face&            theFace,
-                                         std::list<TopoDS_Edge>& theEdges,
+                                         const std::list<TopoDS_Edge>& theEdges,
                                          SMESH_Mesh*                   theMesh,
                                          const bool                    theIsForward,
                                          const bool                    theIgnoreMediumNodes,
@@ -125,7 +125,7 @@ StdMeshers_FaceSide::StdMeshers_FaceSide(const TopoDS_Face&      theFace,
   SMESHDS_Mesh* meshDS = myProxyMesh->GetMeshDS();
 
   int nbDegen = 0;
-  std::list<TopoDS_Edge>::iterator edge = theEdges.begin();
+  std::list<TopoDS_Edge>::const_iterator edge = theEdges.begin();
   for ( int index = 0; edge != theEdges.end(); ++index, ++edge )
   {
     int i = theIsForward ? index : nbEdges-index-1;

src/StdMeshers/StdMeshers_FaceSide.hxx

@@ -83,7 +83,7 @@ public:
    * \brief Wrap several edges. Edges must be properly ordered and oriented.
    */
   StdMeshers_FaceSide(const TopoDS_Face&            theFace,
-                      std::list<TopoDS_Edge>& theEdges,
+                      const std::list<TopoDS_Edge>& theEdges,
                       SMESH_Mesh*                   theMesh,
                       const bool                    theIsForward,
                       const bool                    theIgnoreMediumNodes,
@@ -121,7 +121,7 @@ public:
       ( new StdMeshers_FaceSide( Face,Edge,Mesh,IsForward,IgnoreMediumNodes,FaceHelper,ProxyMesh ));
   }
   static StdMeshers_FaceSidePtr New (const TopoDS_Face&            Face,
-                                     std::list<TopoDS_Edge>& Edges,
+                                     const std::list<TopoDS_Edge>& Edges,
                                      SMESH_Mesh*                   Mesh,
                                      const bool                    IsForward,
                                      const bool                    IgnoreMediumNodes,
@@ -141,9 +141,11 @@ public:
       ( new StdMeshers_FaceSide( Side,Node,Pnt2d1,Pnt2d2,C2d,UFirst,ULast ));
   }
   static StdMeshers_FaceSidePtr New (UVPtStructVec&     theSideNodes,
-                                     const TopoDS_Face& theFace = TopoDS_Face())
+                                     const TopoDS_Face& theFace = TopoDS_Face(),
+                                     const TopoDS_Edge& theEdge = TopoDS_Edge(),
+                                     SMESH_Mesh*        theMesh = 0)
   {
-    return StdMeshers_FaceSidePtr( new StdMeshers_FaceSide( theSideNodes, theFace ));
+    return StdMeshers_FaceSidePtr( new StdMeshers_FaceSide( theSideNodes, theFace, theEdge, theMesh ));
   }
 
   /*!

src/StdMeshers/StdMeshers_Prism_3D.cxx

@@ -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>
@@ -1173,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;
@@ -1181,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 );
+    myHelper->SetElementsOnShape( true );
-    allowHighBndError = !isSimpleBottom( thePrism );
+
+    // 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 );
   }
 
-  if ( !myUseBlock && sweeper.ComputeNodes( *myHelper, tol, allowHighBndError ))
+  if ( myUseBlock ) // use block approach
-  {
-  }
-  else // use block approach
   {
     // loop on nodes inside the bottom face
     Prism_3D::TNode prevBNode;
@@ -2400,7 +2423,7 @@ double StdMeshers_Prism_3D::getSweepTolerance( const Prism_3D::TPrismTopo& thePr
 
 bool StdMeshers_Prism_3D::isSimpleBottom( const Prism_3D::TPrismTopo& thePrism )
 {
-  if ( thePrism.myBottomEdges.size() != 4 )
+  if ( thePrism.myNbEdgesInWires.front() != 4 )
     return false;
 
   // analyse angles between edges
@@ -2432,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)
@@ -4732,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
  */
 //================================================================================
 
@@ -5012,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
+  }
+}

src/StdMeshers/StdMeshers_Prism_3D.hxx

@@ -29,15 +29,12 @@
 
 #include "SMESH_StdMeshers.hxx"
 
-#include "SMDS_MeshNode.hxx"
-#include "SMDS_TypeOfPosition.hxx"
 #include "SMESHDS_Mesh.hxx"
-#include "SMESH_Algo.hxx"
 #include "SMESH_Block.hxx"
-#include "SMESH_Comment.hxx"
 #include "SMESH_Mesh.hxx"
 #include "SMESH_MesherHelper.hxx"
 #include "SMESH_TypeDefs.hxx"
+#include "SMESH_Comment.hxx"
 #include "SMESH_subMesh.hxx"
 
 #include <Adaptor2d_Curve2d.hxx>
@@ -421,10 +418,21 @@ struct StdMeshers_Sweeper
   std::vector< TNodeColumn* > myBndColumns; // boundary nodes
   std::vector< TNodeColumn* > myIntColumns; // internal nodes
 
+  typedef std::vector< double > TZColumn;
+  std::vector< TZColumn > myZColumns; // Z distribution of boundary nodes
+
   bool ComputeNodes( SMESH_MesherHelper& helper,
                      const double        tol,
                      const bool          allowHighBndError );
 
+  bool CheckSameZ();
+
+  bool ComputeNodesOnStraightSameZ( SMESH_MesherHelper& helper );
+
+  bool ComputeNodesOnStraight( SMESH_MesherHelper& helper,
+                               const TopoDS_Face&  bottom,
+                               const TopoDS_Face&  top);
+
 private:
 
   gp_XYZ bndPoint( int iP, int z ) const
@@ -446,6 +454,9 @@ private:
                                  const double                 r,
                                  std::vector< gp_XYZ >&       toIntPoints,
                                  std::vector< double >&       int2BndDist);
+
+  static void fillZColumn( TZColumn&    zColumn,
+                           TNodeColumn& nodes );
 };
 
 // ===============================================
@@ -545,6 +556,11 @@ public:
    */
   bool isSimpleBottom( const Prism_3D::TPrismTopo& thePrism );
 
+  /*!
+   * \brief Defines if all "vertical" EDGEs are straight
+   */
+  bool allVerticalEdgesStraight( const Prism_3D::TPrismTopo& thePrism );
+
   /*!
    * \brief Project mesh faces from a source FACE of one prism to
    *        a source FACE of another prism

src/StdMeshers/StdMeshers_ProjectionUtils.cxx

@@ -2796,9 +2796,7 @@ namespace StdMeshers_ProjectionUtils
    */
   //================================================================================
 
-  typedef list< pair< const SMDS_MeshNode*, const BRepMesh_Triangle* > > TNodeTriaList;
+  void Morph::AddCloseNodes( const SMDS_MeshNode*     srcNode,
-
-  void addCloseNodes( const SMDS_MeshNode*     srcNode,
                              const BRepMesh_Triangle* bmTria,
                              const int                srcFaceID,
                              TNodeTriaList &          noTriQueue )
@@ -2813,7 +2811,7 @@ namespace StdMeshers_ProjectionUtils
         for ( int i = 0, nb = elem->NbNodes(); i < nb; ++i )
         {
           const SMDS_MeshNode* n = elem->GetNode( i );
-          if ( !n->isMarked() )
+          if ( !n->isMarked() /*&& n->getshapeId() == srcFaceID*/ )
             noTriQueue.push_back( make_pair( n, bmTria ));
         }
       }
@@ -2827,10 +2825,10 @@ namespace StdMeshers_ProjectionUtils
    */
   //================================================================================
 
-  const BRepMesh_Triangle* findTriangle( const gp_XY&                            uv,
+  const BRepMesh_Triangle* Morph::FindTriangle( const gp_XY&             uv,
                                                 const BRepMesh_Triangle* bmTria,
-                                         Handle(BRepMesh_DataStructureOfDelaun)& triaDS,
+                                                double                   bc[3],
-                                         double                                  bc[3] )
+                                                int                      triaNodes[3] )
   {
     int   nodeIDs[3];
     gp_XY nodeUVs[3];
@@ -2841,10 +2839,10 @@ namespace StdMeshers_ProjectionUtils
     {
       // check bmTria
 
-      triaDS->ElementNodes( *bmTria, nodeIDs );
+      _triaDS->ElementNodes( *bmTria, nodeIDs );
-      nodeUVs[0] = triaDS->GetNode( nodeIDs[0] ).Coord();
+      nodeUVs[0] = _triaDS->GetNode( nodeIDs[0] ).Coord();
-      nodeUVs[1] = triaDS->GetNode( nodeIDs[1] ).Coord();
+      nodeUVs[1] = _triaDS->GetNode( nodeIDs[1] ).Coord();
-      nodeUVs[2] = triaDS->GetNode( nodeIDs[2] ).Coord();
+      nodeUVs[2] = _triaDS->GetNode( nodeIDs[2] ).Coord();
 
       SMESH_MeshAlgos::GetBarycentricCoords( uv,
                                              nodeUVs[0], nodeUVs[1], nodeUVs[2],
@@ -2852,6 +2850,9 @@ namespace StdMeshers_ProjectionUtils
       if ( bc[0] >= 0 && bc[1] >= 0 && bc[0] + bc[1] <= 1 )
       {
         bc[2] = 1 - bc[0] - bc[1];
+        triaNodes[0] = nodeIDs[0];
+        triaNodes[1] = nodeIDs[1];
+        triaNodes[2] = nodeIDs[2];
         return bmTria;
       }
 
@@ -2862,19 +2863,19 @@ namespace StdMeshers_ProjectionUtils
       gp_XY seg = uv - gc;
 
       bmTria->Edges( linkIDs, ori );
-      int triaID = triaDS->IndexOf( *bmTria );
+      int triaID = _triaDS->IndexOf( *bmTria );
       bmTria = 0;
 
       for ( int i = 0; i < 3; ++i )
       {
-        const BRepMesh_PairOfIndex & triIDs = triaDS->ElementsConnectedTo( linkIDs[i] );
+        const BRepMesh_PairOfIndex & triIDs = _triaDS->ElementsConnectedTo( linkIDs[i] );
         if ( triIDs.Extent() < 2 )
           continue; // no neighbor triangle
 
         // check if a link intersects gc2uv
-        const BRepMesh_Edge & link = triaDS->GetLink( linkIDs[i] );
+        const BRepMesh_Edge & link = _triaDS->GetLink( linkIDs[i] );
-        const BRepMesh_Vertex & n1 = triaDS->GetNode( link.FirstNode() );
+        const BRepMesh_Vertex & n1 = _triaDS->GetNode( link.FirstNode() );
-        const BRepMesh_Vertex & n2 = triaDS->GetNode( link.LastNode() );
+        const BRepMesh_Vertex & n2 = _triaDS->GetNode( link.LastNode() );
         gp_XY uv1 = n1.Coord();
         gp_XY lin = n2.Coord() - uv1; // link direction
 
@@ -2885,7 +2886,7 @@ namespace StdMeshers_ProjectionUtils
         double uSeg = ( uv1 - gc ) ^ lin / crossSegLin;
         if ( 0. <= uSeg && uSeg <= 1. )
         {
-          bmTria = & triaDS->GetElement( triIDs.Index( 1 + ( triIDs.Index(1) == triaID )));
+          bmTria = & _triaDS->GetElement( triIDs.Index( 1 + ( triIDs.Index(1) == triaID )));
           break;
         }
       }
@@ -2893,6 +2894,22 @@ namespace StdMeshers_ProjectionUtils
     return bmTria;
   }
 
+  //================================================================================
+  /*!
+   * \brief Return a triangle sharing a given boundary node
+   *  \param [in] iBndNode - index of the boundary node
+   *  \return const BRepMesh_Triangle* - a found triangle
+   */
+  //================================================================================
+
+  const BRepMesh_Triangle* Morph::GetTriangleNear( int iBndNode )
+  {
+    const BRepMesh::ListOfInteger & linkIds = _triaDS->LinksConnectedTo( iBndNode );
+    const BRepMesh_PairOfIndex &    triaIds = _triaDS->ElementsConnectedTo( linkIds.First() );
+    const BRepMesh_Triangle&           tria = _triaDS->GetElement( triaIds.Index(1) );
+    return &tria;
+  }
+
   //================================================================================
   /*!
    * \brief triangulate the srcFace in 2D
@@ -2997,19 +3014,6 @@ namespace StdMeshers_ProjectionUtils
     const SMDS_MeshNode *srcNode, *tgtNode;
     const BRepMesh_Triangle *bmTria;
 
-    // initialize a queue of nodes with starting triangles
-    TNodeTriaList noTriQueue;
-    size_t iBndSrcN = 1;
-    for ( ; iBndSrcN < _bndSrcNodes.size() &&  noTriQueue.empty();  ++iBndSrcN )
-    {
-      // get a triangle
-      const BRepMesh::ListOfInteger & linkIds = _triaDS->LinksConnectedTo( iBndSrcN );
-      const BRepMesh_PairOfIndex &    triaIds = _triaDS->ElementsConnectedTo( linkIds.First() );
-      const BRepMesh_Triangle&           tria = _triaDS->GetElement( triaIds.Index(1) );
-
-      addCloseNodes( _bndSrcNodes[ iBndSrcN ], &tria, srcFaceID, noTriQueue );
-    }
-
     // un-mark internal src nodes; later we will mark moved nodes
     int nbSrcNodes = 0;
     SMDS_NodeIteratorPtr nIt = _srcSubMesh->GetSubMeshDS()->GetNodes();
@@ -3033,6 +3037,10 @@ namespace StdMeshers_ProjectionUtils
     bool   checkUV = true;
     const SMDS_FacePosition* pos;
 
+    // a queue of nodes with starting triangles
+    TNodeTriaList noTriQueue;
+    size_t iBndSrcN = 1;
+
     while ( nbSrcNodes > 0 )
     {
       while ( !noTriQueue.empty() )
@@ -3048,16 +3056,16 @@ namespace StdMeshers_ProjectionUtils
         // find a delauney triangle containing the src node
         gp_XY uv = tgtHelper.GetNodeUV( srcFace, srcNode, NULL, &checkUV );
         uv *= _scale;
-        bmTria = findTriangle( uv, bmTria, _triaDS, bc );
+        bmTria = FindTriangle( uv, bmTria, bc, nodeIDs );
         if ( !bmTria )
           continue;
 
         // compute new coordinates for a corresponding tgt node
         gp_XY uvNew( 0., 0. ), nodeUV;
-        _triaDS->ElementNodes( *bmTria, nodeIDs );
         for ( int i = 0; i < 3; ++i )
           uvNew += bc[i] * tgtVert( nodeIDs[i]).Coord();
-        uvNew.SetCoord( uvNew.X() / _scale.X(), uvNew.Y() / _scale.Y() );
+        uvNew.SetCoord( uvNew.X() / _scale.X(),
+                        uvNew.Y() / _scale.Y() );
         gp_Pnt xyz = tgtSurface->Value( uvNew );
 
         // find and move tgt node
@@ -3069,7 +3077,7 @@ namespace StdMeshers_ProjectionUtils
         if (( pos = dynamic_cast< const SMDS_FacePosition* >( tgtNode->GetPosition() )))
           const_cast<SMDS_FacePosition*>( pos )->SetParameters( uvNew.X(), uvNew.Y() );
 
-        addCloseNodes( srcNode, bmTria, srcFaceID, noTriQueue );
+        AddCloseNodes( srcNode, bmTria, srcFaceID, noTriQueue );
       }
 
       if ( nbSrcNodes > 0 )
@@ -3077,10 +3085,8 @@ namespace StdMeshers_ProjectionUtils
         // assure that all src nodes are visited
         for ( ; iBndSrcN < _bndSrcNodes.size() &&  noTriQueue.empty();  ++iBndSrcN )
         {
-          const BRepMesh::ListOfInteger & linkIds = _triaDS->LinksConnectedTo( iBndSrcN );
+          const BRepMesh_Triangle* tria = GetTriangleNear( iBndSrcN );
-          const BRepMesh_PairOfIndex &    triaIds = _triaDS->ElementsConnectedTo( linkIds.First() );
+          AddCloseNodes( _bndSrcNodes[ iBndSrcN ], tria, srcFaceID, noTriQueue );
-          const BRepMesh_Triangle&           tria = _triaDS->GetElement( triaIds.Index(1) );
-          addCloseNodes( _bndSrcNodes[ iBndSrcN ], &tria, srcFaceID, noTriQueue );
         }
         if ( noTriQueue.empty() )
         {
@@ -3099,8 +3105,8 @@ namespace StdMeshers_ProjectionUtils
 
   } // Morph::Perform
 
-gp_XY Morph::GetBndUV(const int iNode) const
+  gp_XY Morph::GetBndUV(const int iNode) const
-{
+  {
     return _triaDS->GetNode( iNode ).Coord();
   }
 

src/StdMeshers/StdMeshers_ProjectionUtils.hxx

@@ -159,11 +159,34 @@ namespace StdMeshers_ProjectionUtils
                  const TNodeNodeMap&           src2tgtNodes,
                  const bool                    moveAll);
 
-    // return source boundary nodes. 0-th node is zero
+    // find a triangle containing an UV starting from a given triangle;
+    // return barycentric coordinates of the UV in the found triangle
+    const BRepMesh_Triangle* FindTriangle( const gp_XY&             uv,
+                                           const BRepMesh_Triangle* bmTria,
+                                           double                   bc[3],
+                                           int                      triaNodes[3]);
+
+    // return any delauney triangle neighboring a given boundary node
+    const BRepMesh_Triangle* GetTriangleNear( int iBndNode );
+
+    // return source boundary nodes. 0-th node is NULL so that indices of
+    // boundary nodes correspond to indices used by Delauney mesh
     const std::vector< const SMDS_MeshNode* >& GetBndNodes() const { return _bndSrcNodes; }
 
     // return UV of the i-th source boundary node
     gp_XY GetBndUV(const int iNode) const;
+
+    // return scale factor to convert real UV to/from UV used for Delauney meshing:
+    // delauney_UV = real_UV * scale
+    const gp_XY& GetScale() const { return _scale; }
+
+    typedef std::list< std::pair< const SMDS_MeshNode*, const BRepMesh_Triangle* > > TNodeTriaList;
+
+    // add non-marked nodes surrounding a given one to a list
+    static void AddCloseNodes( const SMDS_MeshNode*     node,
+                               const BRepMesh_Triangle* bmTria,
+                               const int                faceID,
+                               TNodeTriaList &          noTriQueue );
   };
 
   /*!

src/StdMeshers/StdMeshers_Projection_2D.cxx

@@ -1133,7 +1133,7 @@ namespace {
   {
     SMESH_subMesh* faceSM = helper.GetMesh()->GetSubMesh( helper.GetSubShape() );
 
-    if ( helper.IsDistorted2D( faceSM, /*checkUV=*/true ))
+    //if ( helper.IsDistorted2D( faceSM, /*checkUV=*/true ))
     {
       SMESH_MeshEditor editor( helper.GetMesh() );
       SMESHDS_SubMesh* smDS = faceSM->GetSubMeshDS();