Revert "Revert "[EDF31321]: Add bootsrap in smesh's dependence""

This reverts commit 156112ef20db9552d7b91e0a92d295446d13287a.

CMakeLists.txt

@@ -38,9 +38,21 @@ ENDIF()
 
 # Versioning
 # ===========
-SALOME_SETUP_VERSION(9.12.0)
+SALOME_SETUP_VERSION(9.13.0)
 MESSAGE(STATUS "Building ${PROJECT_NAME_UC} ${${PROJECT_NAME_UC}_VERSION} from \"${${PROJECT_NAME_UC}_GIT_SHA1}\"")
 
+
+# Find SalomeBootstrap
+# ===========
+SET(SALOMEBOOTSTRAP_ROOT_DIR $ENV{SALOMEBOOTSTRAP_ROOT_DIR} CACHE PATH "Path to the Salome Bootstrap")
+IF(EXISTS ${SALOMEBOOTSTRAP_ROOT_DIR})
+  FIND_PACKAGE(SalomeBootstrap REQUIRED)
+  ADD_DEFINITIONS(${SALOMEBOOTSTRAP_DEFINITIONS})
+  INCLUDE_DIRECTORIES(${SALOMEBOOTSTRAP_INCLUDE_DIRS})
+ELSE(EXISTS ${SALOMEBOOTSTRAP_ROOT_DIR})
+  MESSAGE(FATAL_ERROR "We absolutely need a Salome Bootstrap, please define SALOMEBOOTSTRAP_ROOT_DIR")
+ENDIF(EXISTS ${SALOMEBOOTSTRAP_ROOT_DIR})
+
 # Find KERNEL
 # ===========
 SET(KERNEL_ROOT_DIR $ENV{KERNEL_ROOT_DIR} CACHE PATH "Path to the Salome KERNEL")

README

@@ -45,7 +45,7 @@ site:
 
 * Check Release Notes and Software Requirements of the latest SALOME release at
 
-  <https://www.salome-platform.org/?page_id=15/>
+  <https://github.com/SalomePlatform/.github/wiki/SAT/>
 
 Note: SALOME Mesh module needs SALOME KERNEL, SALOME GUI and SALOME Geometry as
 pre-requisites.

doc/examples/CMakeLists.txt

@@ -45,13 +45,26 @@ IF(SALOME_BUILD_TESTS)
     SET_TESTS_PROPERTIES(${testname} PROPERTIES ENVIRONMENT "${tests_env}" LABELS "examples")
   ENDFOREACH()
 
+  # Add tests for 'make test'
+  # -------------------------
+  
+  SALOME_GENERATE_TESTS_ENVIRONMENT(tests_env)
+
+  FOREACH(test ${UNSTABLE_TESTS})
+    GET_FILENAME_COMPONENT(testname ${test} NAME_WE)
+    SET(testname "EXAMPLES_${testname}")
+    ADD_TEST(NAME ${testname}
+             COMMAND ${PYTHON_EXECUTABLE} -B ${CMAKE_SOURCE_DIR}/test/test_helper.py ${CMAKE_CURRENT_SOURCE_DIR}/${test})
+    SET_TESTS_PROPERTIES(${testname} PROPERTIES ENVIRONMENT "${tests_env}" LABELS "examples")
+  ENDFOREACH()
+
 ENDIF()
 
 # Install example scripts and data
 # --------------------------------
 
 IF(SALOME_BUILD_TESTS OR SALOME_BUILD_DOC)
-  INSTALL(FILES ${GOOD_TESTS} ${BAD_TESTS} DESTINATION ${EXAMPLES_INSTALL_DIRECTORY})
+  INSTALL(FILES ${GOOD_TESTS} ${BAD_TESTS} ${UNSTABLE_TESTS} DESTINATION ${EXAMPLES_INSTALL_DIRECTORY})
   INSTALL(FILES mechanic.py DESTINATION ${SALOME_INSTALL_BINS})
   INSTALL(DIRECTORY data DESTINATION ${EXAMPLES_INSTALL_DIRECTORY})
 ENDIF()

doc/examples/CTestTestfileInstall.cmake

@@ -32,3 +32,12 @@ FOREACH(tfile ${_all_tests})
   ADD_TEST(${TEST_NAME} python ${PYTHON_TEST_DRIVER} ${TIMEOUT} ${SCRIPTS_DIR}/${tfile})
   SET_TESTS_PROPERTIES(${TEST_NAME} PROPERTIES LABELS "${COMPONENT_NAME};${COMPONENT_NAME}_examples")
 ENDFOREACH()
+
+SET(_unstable_tests ${UNSTABLE_TESTS})
+LIST(SORT _unstable_tests)
+FOREACH(tfile ${_unstable_tests})
+  GET_FILENAME_COMPONENT(BASE_NAME ${tfile} NAME_WE)
+  SET(TEST_NAME UNSTABLE_SMESH_${BASE_NAME})
+  ADD_TEST(${TEST_NAME} python ${PYTHON_TEST_DRIVER} ${TIMEOUT} ${SCRIPTS_DIR}/${tfile})
+  SET_TESTS_PROPERTIES(${TEST_NAME} PROPERTIES LABELS "UNSTABLE;UNSTABLE_SMESH;${COMPONENT_NAME}_examples")
+ENDFOREACH()

doc/examples/tests.set

@@ -152,8 +152,11 @@ ENDIF(NOT WIN32)
 # The following tests can be executed with both 'make test' and 'salome test'.
 # ----------------------------------------------------------------------------
 
-SET(GOOD_TESTS
+SET(UNSTABLE_TESTS
   cartesian_algo.py
+)
+
+SET(GOOD_TESTS
   creating_meshes_ex02.py
   creating_meshes_ex06.py
   creating_meshes_ex07.py

doc/gui/input/1d_meshing_hypo.rst

@@ -1,4 +1,5 @@
 .. _a1d_meshing_hypo_page:
+.. |larr|   unicode:: U+02190 .. LEFTWARDS ARROW
 
 *********************
 1D Meshing Hypotheses
@@ -57,7 +58,7 @@ Adaptive hypothesis
 Arithmetic Progression hypothesis
 #################################
 
-**Arithmetic Progression** hypothesis allows to split edges into segments with a length that changes in arithmetic progression (Lk = Lk-1 + d) beginning from a given starting length and up to a given end length.
+**Arithmetic Progression** hypothesis allows to split edges into segments with a length that changes in arithmetic progression (L\ :sub:`k`\  |larr| L\ :sub:`k-1`\  + d) beginning from a given starting length and up to a given end length.
 
 The splitting direction is defined by the orientation of the underlying geometrical edge. **Reverse Edges** list box allows specifying the edges, for which the splitting should be made in the direction opposite to their orientation. This list box is usable only if a geometry object is selected for meshing. In this case it is possible to select edges to be reversed either directly picking them in the 3D viewer or by selecting the edges or groups of edges in the Object Browser. Use 
 **Add** button to add the selected edges to the list.
@@ -82,7 +83,7 @@ The splitting direction is defined by the orientation of the underlying geometri
 Geometric Progression hypothesis
 ################################
 
-**Geometric Progression** hypothesis allows splitting edges into segments with a length that changes in geometric progression (Lk = Lk-1 * d) starting from a given **Start Length** and with a given **Common Ratio**.
+**Geometric Progression** hypothesis allows splitting edges into segments with a length that changes in geometric progression (L\ :sub:`k`\  |larr| L\ :sub:`k-1`\  * d) starting from a given **Start Length** and with a given **Common Ratio**.
 
 The splitting direction is defined by the orientation of the underlying geometrical edge.
 **Reverse Edges** list box allows specifying the edges, for which the splitting should be made in the direction opposite to their orientation. This list box is usable only if a geometry object is selected for meshing. In this case it is possible to select edges to be reversed either directly picking them in the 3D viewer or by selecting the edges or groups of edges in the Object Browser. Use **Add** button to add the selected edges to the list.
@@ -177,7 +178,7 @@ You can set the type of node distribution for this hypothesis in the **Hypothesi
 
 **Scale Distribution** - length of segments gradually changes depending on the **Scale Factor**, which is a ratio of the first segment length to the last segment length.
 
-Length of segments changes in geometric progression with the common ratio (A) depending on the **Scale Factor** (S) and **Number of Segments** (N) as follows: A = S**(1/(N-1)). For an edge of length L, length of the first segment is L * (1 - A)/(1 - A**N)
+Length of segments changes in geometric progression with the common ratio (A) depending on the **Scale Factor** (S) and **Number of Segments** (N) as follows: A = S\ :sup:`(1/(N-1))`\ . For an edge of length L, length of the first segment is L * (1 - A)/(1 - A\ :sup:`N`\ )
 
 .. image:: ../images/a-nbsegments2.png
 	:align: center

doc/gui/input/constructing_meshes.rst

@@ -210,6 +210,10 @@ When the Preview dialog is closed, the question about the storage of temporarily
 
 These elements can be kept in the mesh.
 
+Since the mesh, its submesh(es) and group(s) are now not empty, their **Object Browser** icons will be updated accordingly.
+
+	.. image:: ../images/image89.jpg
+		:align: center
 
 .. _submesh_order_anchor:
 

doc/gui/input/reload_mesh_from_file.rst

@@ -0,0 +1,41 @@
+.. _reload_mesh_from_file_page:
+
+*******************
+Reload mesh from file
+*******************
+
+This operation allows reload original imported mesh, which was modified by different dedicated tools 
+
+*Reload mesh from file:*
+
+
+#. Select a mesh(es) (and display it in the 3D Viewer if you are going to pick elements by mouse).
+#. From popup menu click on the *Reload from file* item
+
+Each selected mesh will be updated to orignal state with saving display properties
+
+.. image:: ../images/reload_mesh_orig_mesh.png 
+	:align: center
+
+.. centered::
+	Original mesh
+
+|
+
+.. image:: ../images/reload_mesh_modif_mesh.png
+	:align: center
+
+.. centered::
+	Same mesh after apply "Merge Nodes"
+	
+.. image:: ../images/reload_mesh_modif_mesh_props.png
+	:align: center
+
+.. centered::
+	Set display properties
+
+.. image:: ../images/reload_mesh_result.png
+	:align: center
+
+.. centered::
+	Result after reload mesh

idl/SMESH_Gen.idl

@@ -258,6 +258,9 @@ module SMESH
     SMESH_Mesh CreateEmptyMesh()
       raises ( SALOME::SALOME_Exception );
 
+    SMESH_Mesh ReloadMeshFromFile(in SMESH_Mesh sourceMesh)
+      raises(SALOME::SALOME_Exception);
+
     /*!
      * Create Mesh object importing data from given UNV file
      * (UNV supported version is I-DEAS 10)
@@ -302,6 +305,7 @@ module SMESH
                                       out SMESH::DriverMED_ReadStatus theStatus)
       raises (SALOME::SALOME_Exception);
 
+
     /*!
      * Create a dual mesh of a Tetrahedron mesh
      *  \param mesh - TetraHedron mesh to create dual from

src/CMakeLists.txt

@@ -37,6 +37,7 @@ SET(SUBDIRS_COMMON
   SMESHClient
   SMESH_SWIG
   StdMeshers
+  StdMeshers.test
   StdMeshers_I
   SMESH_PY
   Tools

src/Controls/CMakeLists.txt

@@ -51,13 +51,15 @@ SET(_link_LIBRARIES
 # header files / no moc processing
 SET(SMESHControls_HEADERS
   SMESH_ControlsDef.hxx
+  SMESH_ControlsClassifier.hxx
 )
 
 # --- sources ---
 
 # sources / static
 SET(SMESHControls_SOURCES
-  SMESH_Controls.cxx 
+  SMESH_Controls.cxx
+  SMESH_ControlsClassifier.cxx
 )
 
 # --- rules ---

src/Controls/SMESH_Controls.cxx

@@ -40,10 +40,7 @@
 #include <Basics_Utils.hxx>
 
 #include <BRepAdaptor_Surface.hxx>
-#include <BRepBndLib.hxx>
 #include <BRepBuilderAPI_Copy.hxx>
-#include <BRepClass3d_SolidClassifier.hxx>
-#include <BRepClass_FaceClassifier.hxx>
 #include <BRep_Tool.hxx>
 #include <GeomLib_IsPlanarSurface.hxx>
 #include <Geom_CylindricalSurface.hxx>
@@ -1134,6 +1131,7 @@ double AspectRatio3D::GetValue( const TSequenceOfXYZ& P )
     if     (nbNodes==10) nbNodes=4; // quadratic tetrahedron
     else if(nbNodes==13) nbNodes=5; // quadratic pyramid
     else if(nbNodes==15) nbNodes=6; // quadratic pentahedron
+    else if(nbNodes==18) nbNodes=6; // bi-quadratic pentahedron
     else if(nbNodes==20) nbNodes=8; // quadratic hexahedron
     else if(nbNodes==27) nbNodes=8; // tri-quadratic hexahedron
     else return aQuality;
@@ -4493,58 +4491,14 @@ bool ElementsOnSurface::isOnSurface( const SMDS_MeshNode* theNode )
 //  ElementsOnShape
 //================================================================================
 
-namespace {
-  const int theIsCheckedFlag = 0x0000100;
-}
-
-struct ElementsOnShape::Classifier
-{
-  Classifier(): mySolidClfr(0), myProjFace(0), myProjEdge(0), myFlags(0) { myU = myV = 1e100; }
-  ~Classifier();
-  void Init(const TopoDS_Shape& s, double tol, const Bnd_B3d* box = 0 );
-  bool IsOut(const gp_Pnt& p)        { return SetChecked( true ), (this->*myIsOutFun)( p ); }
-  TopAbs_ShapeEnum ShapeType() const { return myShape.ShapeType(); }
-  const TopoDS_Shape& Shape() const  { return myShape; }
-  const Bnd_B3d* GetBndBox() const   { return & myBox; }
-  double Tolerance() const           { return myTol; }
-  bool IsChecked()                   { return myFlags & theIsCheckedFlag; }
-  bool IsSetFlag( int flag ) const   { return myFlags & flag; }
-  void SetChecked( bool is ) { is ? SetFlag( theIsCheckedFlag ) : UnsetFlag( theIsCheckedFlag ); }
-  void SetFlag  ( int flag ) { myFlags |= flag; }
-  void UnsetFlag( int flag ) { myFlags &= ~flag; }
-  void GetParams( double & u, double & v ) const { u = myU; v = myV; }
-
-private:
-  bool isOutOfSolid (const gp_Pnt& p);
-  bool isOutOfBox   (const gp_Pnt& p);
-  bool isOutOfFace  (const gp_Pnt& p);
-  bool isOutOfEdge  (const gp_Pnt& p);
-  bool isOutOfVertex(const gp_Pnt& p);
-  bool isOutOfNone  (const gp_Pnt& /*p*/) { return true; }
-  bool isBox        (const TopoDS_Shape& s);
-
-  TopoDS_Shape prepareSolid( const TopoDS_Shape& theSolid );
-
-  bool (Classifier::*          myIsOutFun)(const gp_Pnt& p);
-  BRepClass3d_SolidClassifier* mySolidClfr;
-  Bnd_B3d                      myBox;
-  GeomAPI_ProjectPointOnSurf*  myProjFace;
-  GeomAPI_ProjectPointOnCurve* myProjEdge;
-  gp_Pnt                       myVertexXYZ;
-  TopoDS_Shape                 myShape;
-  double                       myTol;
-  double                       myU, myV; // result of isOutOfFace() and isOutOfEdge()
-  int                          myFlags;
-};
-
 struct ElementsOnShape::OctreeClassifier : public SMESH_Octree
 {
-  OctreeClassifier( const std::vector< ElementsOnShape::Classifier* >& classifiers );
+  OctreeClassifier( const std::vector< Classifier* >& classifiers );
   OctreeClassifier( const OctreeClassifier*                           otherTree,
-                    const std::vector< ElementsOnShape::Classifier >& clsOther,
-                    std::vector< ElementsOnShape::Classifier >&       cls );
+                    const std::vector< Classifier >& clsOther,
+                    std::vector< Classifier >&       cls );
   void GetClassifiersAtPoint( const gp_XYZ& p,
-                              std::vector< ElementsOnShape::Classifier* >& classifiers );
+                              std::vector< Classifier* >& classifiers );
   size_t GetSize();
 
 protected:
@@ -4553,7 +4507,7 @@ protected:
   void          buildChildrenData();
   Bnd_B3d*      buildRootBox();
 
-  std::vector< ElementsOnShape::Classifier* > myClassifiers;
+  std::vector< Classifier* > myClassifiers;
 };
 
 
@@ -4867,199 +4821,8 @@ bool ElementsOnShape::IsSatisfy (const SMDS_MeshNode* node,
   return !isNodeOut;
 }
 
-void ElementsOnShape::Classifier::Init( const TopoDS_Shape& theShape,
-                                        double              theTol,
-                                        const Bnd_B3d*      theBox )
-{
-  myShape = theShape;
-  myTol   = theTol;
-  myFlags = 0;
-
-  bool isShapeBox = false;
-  switch ( myShape.ShapeType() )
-  {
-  case TopAbs_SOLID:
-  {
-    if (( isShapeBox = isBox( theShape )))
-    {
-      myIsOutFun = & ElementsOnShape::Classifier::isOutOfBox;
-    }
-    else
-    {
-      mySolidClfr = new BRepClass3d_SolidClassifier( prepareSolid( theShape ));
-      myIsOutFun = & ElementsOnShape::Classifier::isOutOfSolid;
-    }
-    break;
-  }
-  case TopAbs_FACE:
-  {
-    Standard_Real u1,u2,v1,v2;
-    Handle(Geom_Surface) surf = BRep_Tool::Surface( TopoDS::Face( theShape ));
-    if ( surf.IsNull() )
-      myIsOutFun = & ElementsOnShape::Classifier::isOutOfNone;
-    else
-    {
-      surf->Bounds( u1,u2,v1,v2 );
-      myProjFace = new GeomAPI_ProjectPointOnSurf;
-      myProjFace->Init( surf, u1,u2, v1,v2, myTol );
-      myIsOutFun = & ElementsOnShape::Classifier::isOutOfFace;
-    }
-    break;
-  }
-  case TopAbs_EDGE:
-  {
-    Standard_Real u1, u2;
-    Handle(Geom_Curve) curve = BRep_Tool::Curve( TopoDS::Edge( theShape ), u1, u2);
-    if ( curve.IsNull() )
-      myIsOutFun = & ElementsOnShape::Classifier::isOutOfNone;
-    else
-    {
-      myProjEdge = new GeomAPI_ProjectPointOnCurve;
-      myProjEdge->Init( curve, u1, u2 );
-      myIsOutFun = & ElementsOnShape::Classifier::isOutOfEdge;
-    }
-    break;
-  }
-  case TopAbs_VERTEX:
-  {
-    myVertexXYZ = BRep_Tool::Pnt( TopoDS::Vertex( theShape ) );
-    myIsOutFun = & ElementsOnShape::Classifier::isOutOfVertex;
-    break;
-  }
-  default:
-    throw SALOME_Exception("Programmer error in usage of ElementsOnShape::Classifier");
-  }
-
-  if ( !isShapeBox )
-  {
-    if ( theBox )
-    {
-      myBox = *theBox;
-    }
-    else
-    {
-      Bnd_Box box;
-      if ( myShape.ShapeType() == TopAbs_FACE )
-      {
-        BRepAdaptor_Surface SA( TopoDS::Face( myShape ), /*useBoundaries=*/false );
-        if ( SA.GetType() == GeomAbs_BSplineSurface )
-          BRepBndLib::AddOptimal( myShape, box,
-                                  /*useTriangulation=*/true, /*useShapeTolerance=*/true );
-      }
-      if ( box.IsVoid() )
-        BRepBndLib::Add( myShape, box );
-      myBox.Clear();
-      myBox.Add( box.CornerMin() );
-      myBox.Add( box.CornerMax() );
-      gp_XYZ halfSize = 0.5 * ( box.CornerMax().XYZ() - box.CornerMin().XYZ() );
-      for ( int iDim = 1; iDim <= 3; ++iDim )
-      {
-        double x = halfSize.Coord( iDim );
-        halfSize.SetCoord( iDim, x + Max( myTol, 1e-2 * x ));
-      }
-      myBox.SetHSize( halfSize );
-    }
-  }
-}
-
-ElementsOnShape::Classifier::~Classifier()
-{
-  delete mySolidClfr; mySolidClfr = 0;
-  delete myProjFace;  myProjFace = 0;
-  delete myProjEdge;  myProjEdge = 0;
-}
-
-TopoDS_Shape ElementsOnShape::Classifier::prepareSolid( const TopoDS_Shape& theSolid )
-{
-  // try to limit tolerance of theSolid down to myTol (issue #19026)
-
-  // check if tolerance of theSolid is more than myTol
-  bool tolIsOk = true; // max tolerance is at VERTEXes
-  for ( TopExp_Explorer exp( theSolid, TopAbs_VERTEX ); exp.More() &&  tolIsOk; exp.Next() )
-    tolIsOk = ( myTol >= BRep_Tool::Tolerance( TopoDS::Vertex( exp.Current() )));
-  if ( tolIsOk )
-    return theSolid;
-
-  // make a copy to prevent the original shape from changes
-  TopoDS_Shape resultShape = BRepBuilderAPI_Copy( theSolid );
-
-  if ( !GEOMUtils::FixShapeTolerance( resultShape, TopAbs_SHAPE, myTol ))
-    return theSolid;
-  return resultShape;
-}
-
-bool ElementsOnShape::Classifier::isOutOfSolid( const gp_Pnt& p )
-{
-  if ( isOutOfBox( p )) return true;
-  mySolidClfr->Perform( p, myTol );
-  return ( mySolidClfr->State() != TopAbs_IN && mySolidClfr->State() != TopAbs_ON );
-}
-
-bool ElementsOnShape::Classifier::isOutOfBox( const gp_Pnt& p )
-{
-  return myBox.IsOut( p.XYZ() );
-}
-
-bool ElementsOnShape::Classifier::isOutOfFace( const gp_Pnt& p )
-{
-  if ( isOutOfBox( p )) return true;
-  myProjFace->Perform( p );
-  if ( myProjFace->IsDone() && myProjFace->LowerDistance() <= myTol )
-  {
-    // check relatively to the face
-    myProjFace->LowerDistanceParameters( myU, myV );
-    gp_Pnt2d aProjPnt( myU, myV );
-    BRepClass_FaceClassifier aClsf ( TopoDS::Face( myShape ), aProjPnt, myTol );
-    if ( aClsf.State() == TopAbs_IN || aClsf.State() == TopAbs_ON )
-      return false;
-  }
-  return true;
-}
-
-bool ElementsOnShape::Classifier::isOutOfEdge( const gp_Pnt& p )
-{
-  if ( isOutOfBox( p )) return true;
-  myProjEdge->Perform( p );
-  bool isOn = ( myProjEdge->NbPoints() > 0 && myProjEdge->LowerDistance() <= myTol );
-  if ( isOn )
-    myU = myProjEdge->LowerDistanceParameter();
-  return !isOn;
-}
-
-bool ElementsOnShape::Classifier::isOutOfVertex( const gp_Pnt& p )
-{
-  return ( myVertexXYZ.Distance( p ) > myTol );
-}
-
-bool ElementsOnShape::Classifier::isBox(const TopoDS_Shape& theShape )
-{
-  TopTools_IndexedMapOfShape vMap;
-  TopExp::MapShapes( theShape, TopAbs_VERTEX, vMap );
-  if ( vMap.Extent() != 8 )
-    return false;
-
-  myBox.Clear();
-  for ( int i = 1; i <= 8; ++i )
-    myBox.Add( BRep_Tool::Pnt( TopoDS::Vertex( vMap( i ))).XYZ() );
-
-  gp_XYZ pMin = myBox.CornerMin(), pMax = myBox.CornerMax();
-  for ( int i = 1; i <= 8; ++i )
-  {
-    gp_Pnt p = BRep_Tool::Pnt( TopoDS::Vertex( vMap( i )));
-    for ( int iC = 1; iC <= 3; ++ iC )
-    {
-      double d1 = Abs( pMin.Coord( iC ) - p.Coord( iC ));
-      double d2 = Abs( pMax.Coord( iC ) - p.Coord( iC ));
-      if ( Min( d1, d2 ) > myTol )
-        return false;
-    }
-  }
-  myBox.Enlarge( myTol );
-  return true;
-}
-
 ElementsOnShape::
-OctreeClassifier::OctreeClassifier( const std::vector< ElementsOnShape::Classifier* >& classifiers )
+OctreeClassifier::OctreeClassifier( const std::vector< Classifier* >& classifiers )
   :SMESH_Octree( new SMESH_TreeLimit )
 {
   myClassifiers = classifiers;
@@ -5068,8 +4831,8 @@ OctreeClassifier::OctreeClassifier( const std::vector< ElementsOnShape::Classifi
 
 ElementsOnShape::
 OctreeClassifier::OctreeClassifier( const OctreeClassifier*                           otherTree,
-                                    const std::vector< ElementsOnShape::Classifier >& clsOther,
-                                    std::vector< ElementsOnShape::Classifier >&       cls )
+                                    const std::vector< Classifier >& clsOther,
+                                    std::vector< Classifier >&       cls )
   :SMESH_Octree( new SMESH_TreeLimit )
 {
   myBox = new Bnd_B3d( *otherTree->getBox() );
@@ -5095,7 +4858,7 @@ OctreeClassifier::OctreeClassifier( const OctreeClassifier*
 
 void ElementsOnShape::
 OctreeClassifier::GetClassifiersAtPoint( const gp_XYZ& point,
-                                         std::vector< ElementsOnShape::Classifier* >& result )
+                                         std::vector< Classifier* >& result )
 {
   if ( getBox()->IsOut( point ))
     return;

src/Controls/SMESH_ControlsClassifier.cxx

@@ -0,0 +1,233 @@
+// Copyright (C) 2007-2024  CEA, EDF, OPEN CASCADE
+//
+// Copyright (C) 2003-2007  OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
+// CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS
+//
+// This library is free software; you can redistribute it and/or
+// modify it under the terms of the GNU Lesser General Public
+// License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+//
+// This library is distributed in the hope that it will be useful,
+// but WITHOUT ANY WARRANTY; without even the implied warranty of
+// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+// Lesser General Public License for more details.
+//
+// You should have received a copy of the GNU Lesser General Public
+// License along with this library; if not, write to the Free Software
+// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
+//
+// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
+//
+
+#include "SMESH_ControlsClassifier.hxx"
+
+#include <GEOMUtils.hxx>
+#include <Utils_SALOME_Exception.hxx>
+
+#include <BRepClass3d_SolidClassifier.hxx>
+#include <BRepClass_FaceClassifier.hxx>
+#include <BRepBuilderAPI_Copy.hxx>
+#include <BRep_Tool.hxx>
+#include <BRepBndLib.hxx>
+#include <TopoDS.hxx>
+#include <TopExp.hxx>
+#include <GeomAPI_ProjectPointOnSurf.hxx>
+#include <GeomAPI_ProjectPointOnCurve.hxx>
+
+namespace SMESH
+{
+    namespace Controls
+    {
+        void Classifier::Init( const TopoDS_Shape& theShape,
+                                                double              theTol,
+                                                const Bnd_B3d*      theBox )
+        {
+          myShape = theShape;
+          myTol   = theTol;
+          myFlags = 0;
+
+          bool isShapeBox = false;
+          switch ( myShape.ShapeType() )
+          {
+          case TopAbs_SOLID:
+          {
+            if (( isShapeBox = isBox( theShape )))
+            {
+              myIsOutFun = & Classifier::isOutOfBox;
+            }
+            else
+            {
+              mySolidClfr = new BRepClass3d_SolidClassifier( prepareSolid( theShape ));
+              myIsOutFun = & Classifier::isOutOfSolid;
+            }
+            break;
+          }
+          case TopAbs_FACE:
+          {
+            Standard_Real u1,u2,v1,v2;
+            Handle(Geom_Surface) surf = BRep_Tool::Surface( TopoDS::Face( theShape ));
+            if ( surf.IsNull() )
+              myIsOutFun = & Classifier::isOutOfNone;
+            else
+            {
+              surf->Bounds( u1,u2,v1,v2 );
+              myProjFace = new GeomAPI_ProjectPointOnSurf;
+              myProjFace->Init( surf, u1,u2, v1,v2, myTol );
+              myIsOutFun = & Classifier::isOutOfFace;
+            }
+            break;
+          }
+          case TopAbs_EDGE:
+          {
+            Standard_Real u1, u2;
+            Handle(Geom_Curve) curve = BRep_Tool::Curve( TopoDS::Edge( theShape ), u1, u2);
+            if ( curve.IsNull() )
+              myIsOutFun = & Classifier::isOutOfNone;
+            else
+            {
+              myProjEdge = new GeomAPI_ProjectPointOnCurve;
+              myProjEdge->Init( curve, u1, u2 );
+              myIsOutFun = & Classifier::isOutOfEdge;
+            }
+            break;
+          }
+          case TopAbs_VERTEX:
+          {
+            myVertexXYZ = BRep_Tool::Pnt( TopoDS::Vertex( theShape ) );
+            myIsOutFun = & Classifier::isOutOfVertex;
+            break;
+          }
+          default:
+            throw SALOME_Exception("Programmer error in usage of Classifier");
+          }
+
+          if ( !isShapeBox )
+          {
+            if ( theBox )
+            {
+              myBox = *theBox;
+            }
+            else
+            {
+              Bnd_Box box;
+              if ( myShape.ShapeType() == TopAbs_FACE )
+              {
+                BRepAdaptor_Surface SA( TopoDS::Face( myShape ), /*useBoundaries=*/false );
+                if ( SA.GetType() == GeomAbs_BSplineSurface )
+                  BRepBndLib::AddOptimal( myShape, box,
+                                          /*useTriangulation=*/true, /*useShapeTolerance=*/true );
+              }
+              if ( box.IsVoid() )
+                BRepBndLib::Add( myShape, box );
+              myBox.Clear();
+              myBox.Add( box.CornerMin() );
+              myBox.Add( box.CornerMax() );
+              gp_XYZ halfSize = 0.5 * ( box.CornerMax().XYZ() - box.CornerMin().XYZ() );
+              for ( int iDim = 1; iDim <= 3; ++iDim )
+              {
+                double x = halfSize.Coord( iDim );
+                halfSize.SetCoord( iDim, x + Max( myTol, 1e-2 * x ));
+              }
+              myBox.SetHSize( halfSize );
+            }
+          }
+        }
+
+        Classifier::~Classifier()
+        {
+          delete mySolidClfr; mySolidClfr = 0;
+          delete myProjFace;  myProjFace = 0;
+          delete myProjEdge;  myProjEdge = 0;
+        }
+
+        TopoDS_Shape Classifier::prepareSolid( const TopoDS_Shape& theSolid )
+        {
+          // try to limit tolerance of theSolid down to myTol (issue #19026)
+
+          // check if tolerance of theSolid is more than myTol
+          bool tolIsOk = true; // max tolerance is at VERTEXes
+          for ( TopExp_Explorer exp( theSolid, TopAbs_VERTEX ); exp.More() &&  tolIsOk; exp.Next() )
+            tolIsOk = ( myTol >= BRep_Tool::Tolerance( TopoDS::Vertex( exp.Current() )));
+          if ( tolIsOk )
+            return theSolid;
+
+          // make a copy to prevent the original shape from changes
+          TopoDS_Shape resultShape = BRepBuilderAPI_Copy( theSolid );
+
+          if ( !GEOMUtils::FixShapeTolerance( resultShape, TopAbs_SHAPE, myTol ))
+            return theSolid;
+          return resultShape;
+        }
+
+        bool Classifier::isOutOfSolid( const gp_Pnt& p )
+        {
+          if ( isOutOfBox( p )) return true;
+          mySolidClfr->Perform( p, myTol );
+          return ( mySolidClfr->State() != TopAbs_IN && mySolidClfr->State() != TopAbs_ON );
+        }
+
+        bool Classifier::isOutOfBox( const gp_Pnt& p )
+        {
+          return myBox.IsOut( p.XYZ() );
+        }
+
+        bool Classifier::isOutOfFace( const gp_Pnt& p )
+        {
+          if ( isOutOfBox( p )) return true;
+          myProjFace->Perform( p );
+          if ( myProjFace->IsDone() && myProjFace->LowerDistance() <= myTol )
+          {
+            // check relatively to the face
+            myProjFace->LowerDistanceParameters( myU, myV );
+            gp_Pnt2d aProjPnt( myU, myV );
+            BRepClass_FaceClassifier aClsf ( TopoDS::Face( myShape ), aProjPnt, myTol );
+            if ( aClsf.State() == TopAbs_IN || aClsf.State() == TopAbs_ON )
+              return false;
+          }
+          return true;
+        }
+
+        bool Classifier::isOutOfEdge( const gp_Pnt& p )
+        {
+          if ( isOutOfBox( p )) return true;
+          myProjEdge->Perform( p );
+          bool isOn = ( myProjEdge->NbPoints() > 0 && myProjEdge->LowerDistance() <= myTol );
+          if ( isOn )
+            myU = myProjEdge->LowerDistanceParameter();
+          return !isOn;
+        }
+
+        bool Classifier::isOutOfVertex( const gp_Pnt& p )
+        {
+          return ( myVertexXYZ.Distance( p ) > myTol );
+        }
+
+        bool Classifier::isBox(const TopoDS_Shape& theShape )
+        {
+          TopTools_IndexedMapOfShape vMap;
+          TopExp::MapShapes( theShape, TopAbs_VERTEX, vMap );
+          if ( vMap.Extent() != 8 )
+            return false;
+
+          myBox.Clear();
+          for ( int i = 1; i <= 8; ++i )
+            myBox.Add( BRep_Tool::Pnt( TopoDS::Vertex( vMap( i ))).XYZ() );
+
+          gp_XYZ pMin = myBox.CornerMin(), pMax = myBox.CornerMax();
+          for ( int i = 1; i <= 8; ++i )
+          {
+            gp_Pnt p = BRep_Tool::Pnt( TopoDS::Vertex( vMap( i )));
+            for ( int iC = 1; iC <= 3; ++ iC )
+            {
+              double d1 = Abs( pMin.Coord( iC ) - p.Coord( iC ));
+              double d2 = Abs( pMax.Coord( iC ) - p.Coord( iC ));
+              if ( Min( d1, d2 ) > myTol )
+                return false;
+            }
+          }
+          myBox.Enlarge( myTol );
+          return true;
+        }
+    }
+}

src/Controls/SMESH_ControlsClassifier.hxx

@@ -0,0 +1,83 @@
+// Copyright (C) 2007-2024  CEA, EDF, OPEN CASCADE
+//
+// Copyright (C) 2003-2007  OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
+// CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS
+//
+// This library is free software; you can redistribute it and/or
+// modify it under the terms of the GNU Lesser General Public
+// License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+//
+// This library is distributed in the hope that it will be useful,
+// but WITHOUT ANY WARRANTY; without even the implied warranty of
+// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+// Lesser General Public License for more details.
+//
+// You should have received a copy of the GNU Lesser General Public
+// License along with this library; if not, write to the Free Software
+// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
+//
+// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
+//
+
+#ifndef _SMESH_CONTROLSCLASSIFIER_HXX_
+#define _SMESH_CONTROLSCLASSIFIER_HXX_
+
+#include "SMESH_Controls.hxx"
+
+#include <TopoDS_Shape.hxx>
+#include <Bnd_B3d.hxx>
+
+class BRepClass3d_SolidClassifier;
+class GeomAPI_ProjectPointOnSurf;
+class GeomAPI_ProjectPointOnCurve;
+
+namespace SMESH
+{
+    namespace Controls
+    {
+        struct SMESHCONTROLS_EXPORT Classifier
+        {
+            Classifier(): mySolidClfr(0), myProjFace(0), myProjEdge(0), myFlags(0) { myU = myV = 1e100; }
+            ~Classifier();
+            void Init(const TopoDS_Shape& s, double tol, const Bnd_B3d* box = 0 );
+            bool IsOut(const gp_Pnt& p)        { return SetChecked( true ), (this->*myIsOutFun)( p ); }
+            TopAbs_ShapeEnum ShapeType() const { return myShape.ShapeType(); }
+            const TopoDS_Shape& Shape() const  { return myShape; }
+            const Bnd_B3d* GetBndBox() const   { return & myBox; }
+            double Tolerance() const           { return myTol; }
+            bool IsChecked()                   { return myFlags & theIsCheckedFlag; }
+            bool IsSetFlag( int flag ) const   { return myFlags & flag; }
+            void SetChecked( bool is ) { is ? SetFlag( theIsCheckedFlag ) : UnsetFlag( theIsCheckedFlag ); }
+            void SetFlag  ( int flag ) { myFlags |= flag; }
+            void UnsetFlag( int flag ) { myFlags &= ~flag; }
+            void GetParams( double & u, double & v ) const { u = myU; v = myV; }
+
+            private:
+            bool isOutOfSolid (const gp_Pnt& p);
+            bool isOutOfBox   (const gp_Pnt& p);
+            bool isOutOfFace  (const gp_Pnt& p);
+            bool isOutOfEdge  (const gp_Pnt& p);
+            bool isOutOfVertex(const gp_Pnt& p);
+            bool isOutOfNone  (const gp_Pnt& /*p*/) { return true; }
+            bool isBox        (const TopoDS_Shape& s);
+
+            TopoDS_Shape prepareSolid( const TopoDS_Shape& theSolid );
+
+            bool (Classifier::*          myIsOutFun)(const gp_Pnt& p);
+            BRepClass3d_SolidClassifier* mySolidClfr;
+            Bnd_B3d                      myBox;
+            GeomAPI_ProjectPointOnSurf*  myProjFace;
+            GeomAPI_ProjectPointOnCurve* myProjEdge;
+            gp_Pnt                       myVertexXYZ;
+            TopoDS_Shape                 myShape;
+            double                       myTol;
+            double                       myU, myV; // result of isOutOfFace() and isOutOfEdge()
+            int                          myFlags;
+
+            static const int theIsCheckedFlag = 0x0000100;
+        };
+    }
+}
+
+#endif

src/Controls/SMESH_ControlsDef.hxx

@@ -26,6 +26,7 @@
 #include "SMESH_Controls.hxx"
 
 #include "SMESH_TypeDefs.hxx"
+#include "SMESH_ControlsClassifier.hxx"
 
 #include <Bnd_B3d.hxx>
 #include <GeomAPI_ProjectPointOnCurve.hxx>
@@ -988,7 +989,6 @@ namespace SMESH{
 
     private:
 
-      struct Classifier;
       struct OctreeClassifier;
 
       void clearClassifiers();

src/DriverCGNS/DriverStructuredCGNS_Write.cxx

@@ -270,17 +270,34 @@ Driver_Mesh::Status DriverStructuredCGNS_Write::Perform()
           if(pRange[4] < pRange[1]) {std::swap(pRange[1],pRange[4]);}
           if(pRange[5] < pRange[2]) {std::swap(pRange[2],pRange[5]);}          
           int cgIndexBoco = 0;              
-          std::string boundaryName = boundaryNames[faceId].empty() ? zoneName + "_" + std::to_string(faceId) : boundaryNames[faceId];
+          std::string boundaryName = boundaryNames[faceId];
+          std::string boundaryNameOrig(boundaryName);
+          bool isInGroup = false;
+          if (boundaryName.empty())
+            boundaryName = zoneName + "_" + std::to_string(faceId);
+          else
+            isInGroup = true;
           if ( bNames.count(boundaryName)!=0 )
             boundaryName = boundaryName + "_" + std::to_string(faceId);
           bNames.insert( boundaryName );
           if(cg_boco_write(_fn, iBase, iZone, boundaryName.c_str(), CGNS_ENUMV(BCTypeNull), 
                             CGNS_ENUMV(PointRange), 2, &pRange[0], &cgIndexBoco) != CG_OK) return addMessage(cg_get_error(), /*fatal = */true);
+          // write also the group name as FamilyName in the same node
+          if (isInGroup)
+          {
+            // go to the BC node
+            if(cg_goto(_fn, iBase, "Zone_t", iZone, "ZoneBC_t", 1,
+                  "BC_t", cgIndexBoco, "end") != CG_OK)
+              return addMessage(cg_get_error(), /*fatal = */true);
+            // write the family name
+            if(cg_famname_write(boundaryNameOrig.c_str()) != CG_OK)
+              return addMessage(cg_get_error(), /*fatal = */true);
+          }
           faceId++;
         }
         // End write boundary
         
-        // Writte Interfaces
+        // Write Interfaces
         for ( TopExp_Explorer fEx( shape, TopAbs_SOLID ); fEx.More(); fEx.Next() )
         {
           TopoDS_Solid neighbourSolid = TopoDS::Solid(fEx.Current()); 
@@ -382,17 +399,34 @@ Driver_Mesh::Status DriverStructuredCGNS_Write::Perform()
           if(pRange[2] < pRange[0]) {std::swap(pRange[0],pRange[2]);}
           if(pRange[1] < pRange[3]) {std::swap(pRange[1],pRange[3]);}          
           int cgIndexBoco = 0;    
-          std::string boundaryName = boundaryNames[edgeId].empty() ? zoneName + "_" + std::to_string(edgeId) : boundaryNames[edgeId];          
+          std::string boundaryName = boundaryNames[edgeId];
+          std::string boundaryNameOrig(boundaryName);
+          bool isInGroup = false;
+          if (boundaryName.empty())
+            boundaryName = zoneName + "_" + std::to_string(edgeId);
+          else
+            isInGroup = true;
           if ( bNames.count(boundaryName)!=0)
             boundaryName = boundaryName + "_" + std::to_string(edgeId);
           bNames.insert( boundaryName );
           if(cg_boco_write(_fn, iBase, iZone, boundaryName.c_str(), CGNS_ENUMV(BCTypeNull), 
                             CGNS_ENUMV(PointRange), 2, &pRange[0], &cgIndexBoco) != CG_OK) return addMessage(cg_get_error(), /*fatal = */true);
+          // write also the group name as FamilyName in the same node
+          if (isInGroup)
+          {
+            // go to the BC node
+            if(cg_goto(_fn, iBase, "Zone_t", iZone, "ZoneBC_t", 1,
+                  "BC_t", cgIndexBoco, "end") != CG_OK)
+              return addMessage(cg_get_error(), /*fatal = */true);
+            // write the family name
+            if(cg_famname_write(boundaryNameOrig.c_str()) != CG_OK)
+              return addMessage(cg_get_error(), /*fatal = */true);
+          }
           edgeId++;
         }
         // End write Boundary
 
-        // Writte Interfaces
+        // Write Interfaces
         for ( TopExp_Explorer fEx( shape, TopAbs_FACE ); fEx.More(); fEx.Next() )
         {
           TopoDS_Face neighbourFace = TopoDS::Face(fEx.Current());

src/OBJECT/SMESH_DeviceActor.cxx

@@ -43,6 +43,7 @@
 #include <vtkObjectFactory.h>
 #include <vtkShrinkFilter.h>
 #include <vtkShrinkPolyData.h>
+#include <vtkTriangleFilter.h>
 
 #include <vtkProperty.h>
 #include <vtkPolyData.h>
@@ -109,6 +110,7 @@ SMESH_DeviceActor
   myMergeFilter = vtkMergeFilter::New();
 
   myGeomFilter = VTKViewer_GeometryFilter::New();
+  myTriangleFilter = vtkTriangleFilter::New();
 
   myTransformFilter = VTKViewer_TransformFilter::New();
 
@@ -153,6 +155,7 @@ SMESH_DeviceActor
   myFaceOrientation->Delete();
 
   myGeomFilter->Delete();
+  myTriangleFilter->Delete();
 
   myTransformFilter->Delete();
 
@@ -251,8 +254,10 @@ SMESH_DeviceActor
 
     anId++; // 3
     myGeomFilter->SetInputConnection( myPassFilter[ anId ]->GetOutputPort() );
+    myTriangleFilter->SetInputConnection(myGeomFilter->GetOutputPort());
 
     anId++; // 4
+    // myPassFilter[ anId ]->SetInputConnection( myTriangleFilter->GetOutputPort() );
     myPassFilter[ anId ]->SetInputConnection( myGeomFilter->GetOutputPort() );
     myPassFilter[ anId + 1 ]->SetInputConnection( myPassFilter[ anId ]->GetOutputPort() );
 

src/OBJECT/SMESH_DeviceActor.h

@@ -46,6 +46,7 @@ class vtkLookupTable;
 class vtkImplicitBoolean;
 class vtkPassThrough;
 class vtkPlaneCollection;
+class vtkTriangleFilter;
 
 class VTKViewer_Transform;
 class VTKViewer_TransformFilter;
@@ -182,6 +183,7 @@ class SMESHOBJECT_EXPORT SMESH_DeviceActor: public vtkLODActor{
 
   bool myStoreClippingMapping;
   VTKViewer_GeometryFilter *myGeomFilter;
+  vtkTriangleFilter* myTriangleFilter = nullptr;
   VTKViewer_TransformFilter *myTransformFilter;
   std::vector<vtkPassThrough*> myPassFilter;
 

src/SMDS/SMDS_FaceOfNodes.hxx

@@ -32,7 +32,7 @@
 class SMDS_EXPORT SMDS_FaceOfNodes: public SMDS_CellOfNodes
 {
  public:
-  void Print(std::ostream & OS) const;
+  virtual void Print(std::ostream & OS) const override;
   SMDS_FaceOfNodes(const SMDS_MeshNode* node1,
                    const SMDS_MeshNode* node2,
                    const SMDS_MeshNode* node3);

src/SMDS/SMDS_MeshElement.hxx

@@ -142,7 +142,7 @@ public:
 
   SMDS_Mesh* GetMesh() const;
 
-  void Print(std::ostream & OS) const;
+  virtual void Print(std::ostream & OS) const;
 
   friend SMDS_EXPORT std::ostream & operator <<(std::ostream & OS, const SMDS_MeshElement *);
   friend class SMDS_ElementFactory;

src/SMDS/SMDS_MeshNode.hxx

@@ -65,7 +65,7 @@ class SMDS_EXPORT SMDS_MeshNode: public SMDS_MeshElement
   virtual bool IsMediumNode(const SMDS_MeshNode* /*node*/) const  { return false; }
   virtual int  NbCornerNodes() const { return 1; }
 
-  void Print(std::ostream & OS) const;
+  virtual void Print(std::ostream & OS) const override;
 
  private:
 

src/SMDS/SMDS_PolygonalFaceOfNodes.hxx

@@ -47,7 +47,7 @@ class SMDS_EXPORT SMDS_PolygonalFaceOfNodes : public SMDS_CellOfNodes
   virtual int NbEdges() const;
   virtual int NbFaces() const;
 
-  virtual void Print (std::ostream & OS) const;
+  virtual void Print (std::ostream & OS) const override;
 
   virtual const SMDS_MeshNode* GetNode(const int ind) const;
 

src/SMDS/SMDS_VolumeOfNodes.hxx

@@ -62,7 +62,7 @@ class SMDS_EXPORT SMDS_VolumeOfNodes: public SMDS_CellOfNodes
                    const int            nbNodes);
   ~SMDS_VolumeOfNodes();
 
-  void Print(std::ostream & OS) const;
+  virtual void Print(std::ostream & OS) const override;
   int NbFaces() const;
   int NbNodes() const;
   int NbEdges() const;

src/SMDS/SMDS_VolumeTool.cxx

@@ -265,26 +265,26 @@ static int QuadPyram_RE [5][9] = { // REVERSED -> FORWARD (EXTERNAL)
   { 3, 8, 0, 9, 4, 12,3, 4, 4 }}; 
 static int QuadPyram_nbN [] = { 8, 6, 6, 6, 6 };
 
-/*   
-//            + N4
-//           /|\
-//         9/ | \10
-//         /  |  \
-//        /   |   \
-//    N3 +----+----+ N5
-//       |    |11  |
-//       |    |    |
-//       |    +13  |                QUADRATIC
-//       |    |    |                PENTAHEDRON
-//     12+    |    +14
-//       |    |    |
-//       |    |    |
-//       |    + N1 |
-//       |   / \   |               
-//       | 6/   \7 |
-//       | /     \ |
-//       |/       \|
-//    N0 +---------+ N2
+/*
+//            + N4                                 +
+//           /|\                                  /|\
+//         9/ | \10                              + | +
+//         /  |  \                              /  |  \
+//        /   |   \                            /   |   \
+//    N3 +----+----+ N5                       +----+----+
+//       |    |11  |                          |    |    |
+//       |    |    |                          |    |    |   Central nodes
+//       |    +13  |      QUADRATIC           | 16 +    |   of bi-quadratic
+//       |    |    |      PENTAHEDRON         |  + |  + |   PENTAHEDRON
+//     12+    |    +14                        +    | 17 +
+//       |    |    |                          |  18|    |
+//       |    |    |                          |    |    |
+//       |    + N1 |                          |    +    |
+//       |   / \   |                          |   / \   |
+//       | 6/   \7 |                          |  +   +  |
+//       | /     \ |                          | /     \ |
+//       |/       \|                          |/       \|
+//    N0 +---------+ N2                       +---------+
 //            8
 */
 static int QuadPenta_F [5][9] = {  // FORWARD
@@ -301,6 +301,20 @@ static int QuadPenta_RE [5][9] = { // REVERSED -> EXTERNAL
   { 0, 12,3, 11,5, 14,2, 8, 0 }}; 
 static int QuadPenta_nbN [] = { 6, 6, 8, 8, 8 };
 
+static int BiQuadPenta_F[5][9] = {  // FORWARD
+  { 0, 6, 1, 7, 2, 8, 0, 0, 0 },
+  { 3, 11,5, 10,4, 9, 3, 3, 3 },
+  { 0, 12,3, 9, 4, 13,1, 6, 16},   //!
+  { 1, 13,4, 10,5, 14,2, 7, 17},   //!
+  { 0, 8, 2, 14,5, 11,3, 12,18} }; //!
+static int BiQuadPenta_RE[5][9] = { // REVERSED -> EXTERNAL
+  { 0, 8, 2, 7, 1, 6, 0, 0, 0 },
+  { 3, 9, 4, 10,5, 11,3, 3, 3 },
+  { 0, 6, 1, 13,4, 9, 3, 12,17},   //!
+  { 1, 7, 2, 14,5, 10,4, 13,16},   //!
+  { 0, 12,3, 11,5, 14,2, 8, 18} }; //!
+static int BiQuadPenta_nbN[] = { 6, 6, 9, 9, 9 };
+
 /*
 //                 13                                                         
 //         N5+-----+-----+N6                          +-----+-----+
@@ -634,6 +648,14 @@ void SMDS_VolumeTool::Inverse ()
     SWAP_NODES( myVolumeNodes, 9, 11 );
     SWAP_NODES( myVolumeNodes, 13, 14 );
     break;
+  case 18:
+    SWAP_NODES(myVolumeNodes, 1, 2);
+    SWAP_NODES(myVolumeNodes, 4, 5);
+    SWAP_NODES(myVolumeNodes, 6, 8);
+    SWAP_NODES(myVolumeNodes, 9, 11);
+    SWAP_NODES(myVolumeNodes, 13, 14);
+    SWAP_NODES(myVolumeNodes, 16, 17);
+    break;
   case 20:
     SWAP_NODES( myVolumeNodes, 1, 3 );
     SWAP_NODES( myVolumeNodes, 5, 7 );
@@ -677,6 +699,7 @@ SMDS_VolumeTool::VolumeType SMDS_VolumeTool::GetVolumeType() const
   case 10: return QUAD_TETRA;
   case 13: return QUAD_PYRAM;
   case 15: return QUAD_PENTA;
+  case 18: return QUAD_PENTA;
   case 20: return QUAD_HEXA;
   case 27: return QUAD_HEXA;
   default: break;
@@ -862,6 +885,8 @@ double SMDS_VolumeTool::GetSize() const
                            myVolumeNodes[ vtab[i][2] ],
                            myVolumeNodes[ vtab[i][3] ]);
     }
+    if (!myVolForward && V < 0)
+      V *= -1;
   }
   return V;
 }
@@ -1705,6 +1730,16 @@ int SMDS_VolumeTool::GetCenterNodeIndex( int faceIndex ) const
       return faceIndex + 19;
     }
   }
+  else if (myAllFacesNbNodes && myVolumeNodes.size() == 18) // element with 18 nodes
+  {
+    switch (faceIndex) {
+      case 2: return 15;
+      case 3: return 16;
+      case 4: return 17;
+      default:
+        return -2;
+    }
+  }
   return -1;
 }
 
@@ -1727,6 +1762,7 @@ int SMDS_VolumeTool::GetOppFaceIndex( int faceIndex ) const
     switch ( myVolumeNodes.size() ) {
     case 6:
     case 15:
+    case 18:
       if ( faceIndex == 0 || faceIndex == 1 )
         ind = 1 - faceIndex;
       break;
@@ -2469,6 +2505,7 @@ bool SMDS_VolumeTool::setFace( int faceIndex ) const
         myMaxFaceNbNodes         = sizeof(QuadPyram_F[0])/sizeof(QuadPyram_F[0][0]);
         break;
       case 15:
+      case 18:
         myAllFacesNodeIndices_F  = &QuadPenta_F [0][0];
         //myAllFacesNodeIndices_FE = &QuadPenta_FE[0][0];
         myAllFacesNodeIndices_RE = &QuadPenta_RE[0][0];
@@ -2535,7 +2572,8 @@ SMDS_VolumeTool::VolumeType SMDS_VolumeTool::GetType(int nbNodes)
   case 8: return HEXA;
   case 10: return QUAD_TETRA;
   case 13: return QUAD_PYRAM;
-  case 15: return QUAD_PENTA;
+  case 15:
+  case 18: return QUAD_PENTA;
   case 20:
   case 27: return QUAD_HEXA;
   case 12: return HEX_PRISM;

src/SMESH/SMESH_Hypothesis.cxx

@@ -47,8 +47,12 @@ SMESH_Hypothesis::SMESH_Hypothesis(int hypId,
   _type           = PARAM_ALGO;
   _shapeType      = 0;  // to be set by algo with TopAbs_Enum
   _param_algo_dim = -1; // to be set by algo parameter
-  StudyContextStruct* myStudyContext = gen->GetStudyContext();
-  myStudyContext->mapHypothesis[hypId] = this;
+  
+  if ( _gen )
+  {
+    StudyContextStruct* myStudyContext = gen->GetStudyContext();
+    myStudyContext->mapHypothesis[hypId] = this;
+  }    
 }
 
 //=============================================================================
@@ -107,16 +111,18 @@ int SMESH_Hypothesis::GetShapeType() const
 void SMESH_Hypothesis::NotifySubMeshesHypothesisModification()
 {
   // for all meshes in study
-
-  StudyContextStruct* myStudyContext = _gen->GetStudyContext();
-  map<int, SMESH_Mesh*>::iterator itm;
-  for (itm = myStudyContext->mapMesh.begin();
-       itm != myStudyContext->mapMesh.end();
-       itm++)
+  if ( _gen )
   {
-    SMESH_Mesh* mesh = (*itm).second;
-    mesh->NotifySubMeshesHypothesisModification( this );
-  }
+    StudyContextStruct* myStudyContext = _gen->GetStudyContext();
+    map<int, SMESH_Mesh*>::iterator itm;
+    for (itm = myStudyContext->mapMesh.begin();
+        itm != myStudyContext->mapMesh.end();
+        itm++)
+    {
+      SMESH_Mesh* mesh = (*itm).second;
+      mesh->NotifySubMeshesHypothesisModification( this );
+    }
+  }    
 }
 
 //=============================================================================
@@ -148,13 +154,16 @@ void SMESH_Hypothesis::SetLibName(const char* theLibName)
 
 SMESH_Mesh* SMESH_Hypothesis::GetMeshByPersistentID(int id) const
 {
-  StudyContextStruct* myStudyContext = _gen->GetStudyContext();
-  map<int, SMESH_Mesh*>::iterator itm = myStudyContext->mapMesh.begin();
-  for ( ; itm != myStudyContext->mapMesh.end(); itm++)
+  if ( _gen )
   {
-    SMESH_Mesh* mesh = (*itm).second;
-    if ( mesh->GetMeshDS()->GetPersistentId() == id )
-      return mesh;
+    StudyContextStruct* myStudyContext = _gen->GetStudyContext();
+    map<int, SMESH_Mesh*>::iterator itm = myStudyContext->mapMesh.begin();
+    for ( ; itm != myStudyContext->mapMesh.end(); itm++)
+    {
+      SMESH_Mesh* mesh = (*itm).second;
+      if ( mesh->GetMeshDS()->GetPersistentId() == id )
+        return mesh;
+    }
   }
   return 0;
 }

src/SMESH/SMESH_Mesh.cxx

@@ -86,6 +86,8 @@
 namespace fs=boost::filesystem;
 #endif
 
+#include <unordered_set>
+
 // maximum stored group name length in MED file
 #define MAX_MED_GROUP_NAME_LENGTH 80
 
@@ -713,6 +715,8 @@ SMESH_Mesh::AddHypothesis(const TopoDS_Shape & aSubShape,
         }
       }
     }
+
+    ret = CheckHypothesesOnSubMeshes(subMesh, anHyp, event);
   }
   HasModificationsToDiscard(); // to reset _isModified flag if a mesh becomes empty
   GetMeshDS()->Modified();
@@ -1002,6 +1006,78 @@ SMESH_Hypothesis * SMESH_Mesh::GetHypothesis(const int anHypId) const
   return anHyp;
 }
 
+//================================================================================
+/*!
+ * \brief Iterates hypotesis for all sub-meshes of the given sub-mesh and checks
+  algo state with the given event. The goal is to address hypothesis those are
+  not directly affected by changing of an algorithm of the given sub-shape.
+  It is essential to rebuild propagation chains of such hypotheses, otherwise the chains
+  are being cleared after editing of the algorithm and never rebuilt again.
+  * \param subMesh - the main sub-mesh to check sub-meshes of
+  * \param anHyp - the hypothesis changed on the given sub-mesh, we need to skip it from checking
+  * \param event - the given event
+  * \retval SMESH_Hypothesis::Hypothesis_Status - HYP_OK if no errors found, otherwise the most severe error
+ */
+//================================================================================
+SMESH_Hypothesis::Hypothesis_Status SMESH_Mesh::CheckHypothesesOnSubMeshes(
+  SMESH_subMesh* subMesh,
+  const SMESH_Hypothesis* anHyp,
+  const SMESH_subMesh::algo_event event) const
+{
+  SMESH_Hypothesis::Hypothesis_Status ret = SMESH_Hypothesis::Hypothesis_Status::HYP_OK;
+
+  // Cache the processed hypotheses for performance reasons.
+  // Given hypothesis is already processed, so should be skipped.
+  std::unordered_set<const SMESH_Hypothesis*> processedHypotheses = { anHyp };
+
+  // Look through sub-meshes of the given sub-mesh
+  SMESH_subMeshIteratorPtr smIt = subMesh->getDependsOnIterator(false, false);
+  while (smIt->more())
+  {
+    const SMESH_subMesh* sm = smIt->next();
+    const SMESH_Algo* algo = sm->GetAlgo();
+    if (!algo)
+      continue;
+
+    const SMESH_HypoFilter* hypoKind = algo->GetCompatibleHypoFilter(false);
+    if (!hypoKind)
+      continue;
+
+    std::list <const SMESHDS_Hypothesis*> usedHyps;
+    if (!GetHypotheses(sm, *hypoKind, usedHyps, true))
+      continue;
+
+    // Look through hypotheses used by algo
+    for (const auto* usedHyp : usedHyps)
+    {
+      SMESH_Hypothesis* hyp = GetHypothesis(usedHyp->GetID());
+      if (hyp == anHyp)
+        continue;
+
+      if (processedHypotheses.find(hyp) != processedHypotheses.end())
+        continue;
+
+      processedHypotheses.insert(hyp); // Cache the hypothesis pointer
+
+      // Hypoteses restricted by Propagation only because of failed tests.
+      // It's ok for now, because this method was created to fix propagation issue.
+      // It should be investigated more if we find similar issues with other hypotheses.
+      const char* hypName = hyp->GetName();
+      if (strcmp(hypName, "Propagation") != 0)
+        continue;
+
+      const SMESH_Hypothesis::Hypothesis_Status ret2 = subMesh->SubMeshesAlgoStateEngine(event, hyp, true);
+      if (ret2 > ret)
+      {
+        ret = ret2;
+        break;
+      }
+    }
+  }
+
+  return ret;
+}
+
 //=============================================================================
 /*!
  *

src/SMESH/SMESH_Mesh.hxx

@@ -175,6 +175,11 @@ class SMESH_EXPORT SMESH_Mesh
 
   SMESH_Hypothesis * GetHypothesis(const int aHypID) const;
 
+  SMESH_Hypothesis::Hypothesis_Status CheckHypothesesOnSubMeshes(
+    SMESH_subMesh* subMesh,
+    const SMESH_Hypothesis* anHyp,
+    const SMESH_subMesh::algo_event event) const;
+
   const std::list<SMESHDS_Command*> & GetLog();
 
   void ClearLog();

src/SMESH/SMESH_MeshEditor.cxx

@@ -9483,7 +9483,8 @@ void SMESH_MeshEditor::ConvertToQuadratic(const bool theForce3d, const bool theT
       case SMDSEntity_TriQuad_Hexa:
         NewVolume = aHelper.AddVolume(nodes[0], nodes[1], nodes[2], nodes[3],
                                       nodes[4], nodes[5], nodes[6], nodes[7], id, theForce3d);
-        for ( size_t i = 20; i < nodes.size(); ++i ) // rm central nodes
+        for (size_t i = 8; i < nodes.size(); ++i) // rm central nodes from each edge
+        //for (size_t i = 20; i < nodes.size(); ++i) // rm central nodes from each edge
           if ( nodes[i]->NbInverseElements() == 0 )
             GetMeshDS()->RemoveFreeNode( nodes[i], /*sm=*/0, /*fromGroups=*/true );
         break;
@@ -9496,7 +9497,9 @@ void SMESH_MeshEditor::ConvertToQuadratic(const bool theForce3d, const bool theT
       case SMDSEntity_BiQuad_Penta:
         NewVolume = aHelper.AddVolume(nodes[0], nodes[1], nodes[2],
                                       nodes[3], nodes[4], nodes[5], id, theForce3d);
-        for ( size_t i = 15; i < nodes.size(); ++i ) // rm central nodes
+
+        for (size_t i = 6; i < nodes.size(); ++i) // rm central nodes
+        //for ( size_t i = 15; i < nodes.size(); ++i ) // rm central nodes
           if ( nodes[i]->NbInverseElements() == 0 )
             GetMeshDS()->RemoveFreeNode( nodes[i], /*sm=*/0, /*fromGroups=*/true );
         break;
@@ -13037,15 +13040,58 @@ int SMESH_MeshEditor::MakeBoundaryMesh(const TIDSortedElemSet& elements,
           if (iQuad)
             for ( inode = 1; inode < nbFaceNodes; inode += 2)
               nodes.push_back( nn[inode] ); // add medium nodes
-          int iCenter = vTool.GetCenterNodeIndex(iface); // for HEX27
-          if ( iCenter > 0 )
-            nodes.push_back( vTool.GetNodes()[ iCenter ] );
 
-          if (const SMDS_MeshElement * f = aMesh->FindElement( nodes,
-                                                               SMDSAbs_Face, /*noMedium=*/false ))
-            presentBndElems.push_back( f );
+          // for triangle face for Penta18 (BiQuadratic pentahedron) return -2
+          // because we haven't center node on triangle side, but it's need for create biquadratic face
+          int iCenter = vTool.GetCenterNodeIndex(iface); // for HEX27
+
+          // for triangle faces for Penta18 (BiQuadratic pentahedron) firstly check, exist face or not
+          // if not - create node in middle face
+          if (iCenter == -2)
+          {
+            SMDS_ElemIteratorPtr itF = nodes[0]->GetInverseElementIterator(SMDSAbs_Face);
+            bool isFound = false;
+            while (itF->more())
+            {
+              const SMDS_MeshElement* e = itF->next();
+              int nbNodesToCheck = e->NbNodes();
+              if (nbNodesToCheck == (int)nodes.size() + 1)
+              {
+                for (size_t i = 1; e && i < nodes.size() - 1; ++i)
+                {
+                  int nodeIndex = e->GetNodeIndex(nodes[i]);
+                  if (nodeIndex < 0 || nodeIndex >= nbNodesToCheck)
+                    e = 0;
+                }
+                if (e)
+                {
+                  presentBndElems.push_back(e);
+                  isFound = true;
+                }
+              }
+            }
+
+            if (!isFound)
+            {
+              SMESH_MesherHelper aHelper(*myMesh);
+              double bc[3];
+              vTool.GetFaceBaryCenter(iface, bc[0], bc[1], bc[2]);
+              auto aNodeC = aHelper.AddNode(bc[0], bc[1], bc[2]);
+              nodes.push_back(aNodeC);
+              missingBndElems.push_back(nodes);
+            }
+          }
           else
-            missingBndElems.push_back( nodes );
+          {
+            if (iCenter > 0)
+              nodes.push_back(vTool.GetNodes()[iCenter]);
+
+            if (const SMDS_MeshElement* f = aMesh->FindElement(nodes,
+              SMDSAbs_Face, /*noMedium=*/false))
+              presentBndElems.push_back(f);
+            else
+              missingBndElems.push_back(nodes);
+          }
 
           if ( targetMesh != myMesh )
           {

src/SMESHDS/SMESHDS_Script.cxx

@@ -35,7 +35,8 @@ using namespace std;
 //purpose  : 
 //=======================================================================
 SMESHDS_Script::SMESHDS_Script(bool theIsEmbeddedMode):
-  myIsEmbeddedMode(theIsEmbeddedMode)
+  myIsEmbeddedMode(theIsEmbeddedMode),
+  myIsModified(true)
 {
   //cerr << "=========================== myIsEmbeddedMode " << myIsEmbeddedMode << endl;
 }

src/SMESHDS/SMESHDS_SubMesh.hxx

@@ -86,6 +86,8 @@ class SMESHDS_EXPORT SMESHDS_SubMesh : public SMDS_ElementHolder
   virtual void tmpClear();
   virtual void add( const SMDS_MeshElement* element );
   virtual void compact() {}
+  // Commented out to avoid SMESH_netgen_runner_1D2D3D test failure
+  // virtual void clear() override { Clear(); }
 
  private:
 

src/SMESHGUI/SMESHGUI.cxx

@@ -122,6 +122,9 @@
 
 #include "SMESH_version.h"
 
+#include "SMESHDS_Mesh.hxx"
+#include "SMESH_Mesh.hxx"
+
 #include "SMESH_Actor.h"
 #include "SMESH_ActorUtils.h"
 #include "SMESH_Client.hxx"
@@ -221,6 +224,8 @@ namespace
 
   void ExportMeshToFile(int theCommandID);
 
+  void ReloadMeshFromFile(int theCommandID);
+
   void SetDisplayMode(int theCommandID, VTK::MarkerMap& theMarkerMap);
 
   void SetDisplayEntity(int theCommandID);
@@ -410,9 +415,6 @@ namespace
         {
           _PTR(SObject) aMeshSO = SMESH::FindSObject( aMeshes[i] );
           if ( aMeshSO ) {
-            _PTR(StudyBuilder) aBuilder = aStudy->NewBuilder();
-            _PTR(AttributePixMap) aPixmap = aBuilder->FindOrCreateAttribute( aMeshSO, "AttributePixMap" );
-            aPixmap->SetPixMap( "ICON_SMESH_TREE_MESH_IMPORTED" );
             if ( theCommandID == SMESHOp::OpImportUNV ) // mesh names aren't taken from the file for UNV import
               SMESH::SetName( aMeshSO, QFileInfo(filename).fileName() );
 
@@ -607,7 +609,7 @@ namespace
     {
       format = "CGNS";
       notSupportedElemTypes.push_back( SMESH::Entity_Ball );
-    }    
+    }
     else if ( isGMF )
     {
       format = "GMF";
@@ -707,7 +709,7 @@ namespace
       checkBoxes << QObject::tr( "CGNS_EXPORT_ELEMS_BY_TYPE" ) << QObject::tr("STRUCTUREDCGNS");
 
       SalomeApp_CheckFileDlg* fd = new SalomeApp_CheckFileDlg ( SMESHGUI::desktop(), false, checkBoxes, true, true );
-      
+
       fd->setWindowTitle( aTitle );
       fd->setNameFilter( QObject::tr( "CGNS_FILES_FILTER" ) + " (*.cgns)" );
 
@@ -717,18 +719,18 @@ namespace
       SMESHGUI_FileValidator* fv = new SMESHGUI_FileValidator( fd );
       fd->setValidator( fv );
       fd->SetChecked( option, 0 );
-      
+
       if ( fd->exec() )
       {
         aFilename     = fd->selectedFile();
-        structureCGNS = fd->IsChecked( 1 );        
+        structureCGNS = fd->IsChecked( 1 );
       }
 
       toOverwrite = fv->isOverwrite( aFilename );
       option      = fd->IsChecked( 0 );
       SMESHGUI::resourceMgr()->setValue("SMESH", theOptionResource, option );
       toCreateGroups = option;
-      
+
       delete fd;
     }
     else if ( isUNV || isDAT ) // Export to [ UNV | DAT ] - one option
@@ -768,11 +770,11 @@ namespace
       QStringList filters;
       aFilterMap.insert( QObject::tr( "STL_ASCII_FILES_FILTER" ) + " (*.stl)", 1 );
       aFilterMap.insert( QObject::tr( "STL_BIN_FILES_FILTER" )   + " (*.stl)", 0 );
-      
+
       QMap<QString, int>::const_iterator it = aFilterMap.begin();
       for ( ; it != aFilterMap.end(); ++it )
-        filters.push_back( it.key() );    
-     
+        filters.push_back( it.key() );
+
       SUIT_FileDlg* fd = new SUIT_FileDlg( SMESHGUI::desktop(), false, true, true );
       fd->setWindowTitle( aTitle );
       fd->setNameFilters( filters );
@@ -1031,7 +1033,7 @@ namespace
                                     aFilename.toUtf8().data(),
                                     toOverwrite && aMeshIndex == 0,
                                     toCreateGroups );
-            else                                  
+            else
               aMeshItem->ExportStructuredCGNS( aMeshOrGroup,
                                                 aFilename.toUtf8().data(),
                                                 toOverwrite && aMeshIndex == 0 );
@@ -1069,16 +1071,94 @@ namespace
           SUIT_MessageBox::warning(SMESHGUI::desktop(),
                                    QObject::tr("SMESH_WRN_WARNING"),
                                    QObject::tr("SMESH_EXPORT_FAILED") + SalomeApp_Tools::ExceptionToString(S_ex));
-        }          
+        }
         else
           SUIT_MessageBox::warning(SMESHGUI::desktop(),
                                    QObject::tr("SMESH_WRN_WARNING"),
-                                   SalomeApp_Tools::ExceptionToString(S_ex));                                           
+                                   SalomeApp_Tools::ExceptionToString(S_ex));
         wc.resume();
       }
     }
   }
 
+  //================================================================================
+  /*!
+   * \brief Reload selected mesh from file a file
+   */
+   //================================================================================
+  void ReloadMeshFromFile(int theCommandID)
+  {
+    LightApp_SelectionMgr* aSel = SMESHGUI::selectionMgr();
+    SALOME_ListIO selected;
+    if (aSel)
+      aSel->selectedObjects(selected);
+
+    QList< QPair< SMESH::SMESH_IDSource_var, QString > >           aMeshList;
+    QList< QPair< SMESH::SMESH_IDSource_var, QString > >::iterator aMeshIter;
+    SALOME_ListIteratorOfListIO It(selected);
+
+    // Iterate by all selected
+    for (; It.More(); It.Next())
+    {
+      Handle(SALOME_InteractiveObject) anIObject = It.Value();
+      SMESH::SMESH_IDSource_var aMeshItem =
+        SMESH::IObjectToInterface<SMESH::SMESH_IDSource>(anIObject);
+      
+      if (aMeshItem->_is_nil()) {
+        SUIT_MessageBox::warning(SMESHGUI::desktop(),
+          QObject::tr("SMESH_WRN_WARNING"),
+          QObject::tr("SMESH_BAD_MESH_SELECTION"));
+        continue;
+      }
+
+      SMESH::SMESH_Mesh_var aMeshByIO = SMESH::GetMeshByIO(anIObject);
+      SMESH::SelectionProxy aMesh = SMESH::SelectionProxy(aMeshItem);
+      SMESH::MedInfo anInfo = aMesh.medFileInfo();
+      if (!anInfo.isValid())
+        continue;
+
+      SMESH::ListOfGroups aGroups = *aMeshByIO->GetGroups();
+      for (int i = 0; i < aGroups.length(); ++i)
+      {
+        auto X = aGroups[i];
+        _PTR(SObject) aGroupSObject = SMESH::FindSObject(X);
+        SMESH_Actor* anActor = SMESH::FindActorByEntry(aGroupSObject->GetID().c_str());
+        SMESH::smIdType_array aMeshInfo = *X->GetMeshInfo();
+        
+        if (anActor)
+        {
+          vtkTypeBool isVisib = anActor->GetVisibility();
+          SMESH::UpdateView(!isVisib ? SMESH::eDisplay : SMESH::eErase, aGroupSObject->GetID().c_str());
+        }
+        else
+          SMESH::UpdateView(SMESH::eErase, aGroupSObject->GetID().c_str());
+      }
+
+      SMESH::mesh_array_var aMeshes = new SMESH::mesh_array;
+      {
+        // re-import mesh
+        SMESH::SMESH_Mesh_var aReloadedMesh = SMESHGUI::GetSMESHGen()->ReloadMeshFromFile(aMeshByIO);
+
+        QStringList anEntryList;
+
+        _PTR(SObject) aMeshSO = SMESH::FindSObject(aReloadedMesh);
+        if (aMeshSO) {
+          anEntryList.append(aMeshSO->GetID().c_str());
+        }
+        SMESHGUI::GetSMESHGUI()->updateObjBrowser();
+        if (LightApp_Application* anApp =
+          dynamic_cast<LightApp_Application*>(SUIT_Session::session()->activeApplication()))
+          anApp->browseObjects(anEntryList);
+      }
+
+      SMESH::SMESH_GroupBase_var aGroup = SMESH::SMESH_GroupBase::_narrow(aMeshItem);
+      QString aMeshName = anIObject->getName();
+      aMeshList.append(QPair< SMESH::SMESH_IDSource_var, QString >(aMeshItem, aMeshName));
+
+    }
+    SMESH::UpdateView();
+  }
+
   inline void InverseEntityMode(unsigned int& theOutputMode,
                                 unsigned int  theMode)
   {
@@ -2691,7 +2771,12 @@ bool SMESHGUI::OnGUIEvent( int theCommandID )
       ::ImportMeshesFromFile(GetSMESHGen(),theCommandID);
       break;
     }
-
+  case SMESHOp::OpReloadFromFile:
+  {
+    if (isStudyLocked()) break;
+    ::ReloadMeshFromFile(theCommandID);
+    break;
+  }
   case SMESHOp::OpFileInformation:
     {
       SALOME_ListIO selected;
@@ -3484,7 +3569,6 @@ bool SMESHGUI::OnGUIEvent( int theCommandID )
       }
       break;
     }
-
   case SMESHOp::OpFindElementByPoint:
     {
       startOperation( theCommandID );
@@ -4196,7 +4280,7 @@ void SMESHGUI::createSMESHAction( const int id, const QString& po_id, const QStr
     pix = resMgr->loadPixmap( "SMESH", tr( QString( "ICON_%1" ).arg( po_id ).toLatin1().data() ), false );
   if ( !pix.isNull() )
     icon = QIcon( pix );
-  
+
   QString tooltip    = tr( QString( "TOP_%1" ).arg( po_id ).toLatin1().data() ),
           menu       = tr( QString( "MEN_%1" ).arg( po_id ).toLatin1().data() ),
           status_bar = tr( QString( "STB_%1" ).arg( po_id ).toLatin1().data() );
@@ -4320,6 +4404,7 @@ void SMESHGUI::initialize( CAM_Application* app )
   //createSMESHAction( SMESHOp::OpStdInfo, "STD_INFO",        "ICON_STD_INFO" );
   //createSMESHAction( SMESHOp::OpWhatIs, "WHAT_IS",         "ICON_WHAT_IS" ); // VSR: issue #0021242 (eliminate "Mesh Element Information" command)
   createSMESHAction( SMESHOp::OpFindElementByPoint,   "FIND_ELEM",               "ICON_FIND_ELEM" );
+  createSMESHAction( SMESHOp::OpReloadFromFile,       "RELOAD_FROM_FILE");
   //update
   createSMESHAction( SMESHOp::OpFreeNode,              "FREE_NODE",               "ICON_FREE_NODE",     0, true );
   createSMESHAction( SMESHOp::OpEqualNode,             "EQUAL_NODE",              "ICON_EQUAL_NODE",    0, true );
@@ -4932,6 +5017,7 @@ void SMESHGUI::initialize( CAM_Application* app )
   createPopupItem( SMESHOp::OpFileInformation,   OB, mesh, "&& selcount=1 && isImported" );
   createPopupItem( SMESHOp::OpMeshInformation,   OB, mesh_part );
   createPopupItem( SMESHOp::OpFindElementByPoint,OB, mesh_group, "&& selcount=1 && " + hasElems );
+  createPopupItem( SMESHOp::OpReloadFromFile,    OB, mesh, "&& isImported");
   createPopupItem( SMESHOp::OpOverallMeshQuality,OB, mesh_part );
   popupMgr()->insert( separator(), -1, 0 );
   createPopupItem( SMESHOp::OpCreateGroup,       OB, mesh, "&& selcount=1" );

src/SMESHGUI/SMESHGUI_AddMeshElementDlg.cxx

@@ -70,6 +70,8 @@
 #include <vtkPolyDataMapper.h>
 #include <vtkProperty.h>
 #include <vtkCellData.h>
+#include <vtkTriangleFilter.h>
+#include <vtkGeometryFilter.h>
 
 // Qt includes
 #include <QComboBox>
@@ -108,6 +110,8 @@ namespace SMESH
     SALOME_Actor* myFaceOrientation;
     vtkPolyDataMapper* myFaceOrientationDataMapper;
     SMESH_FaceOrientationFilter* myFaceOrientationFilter;
+    vtkSmartPointer<vtkGeometryFilter> myGeometryFilter = nullptr;
+    vtkSmartPointer<vtkTriangleFilter> myTriangleFilter = nullptr;
 
   public:
     TElementSimulation (SalomeApp_Application* theApplication)
@@ -120,9 +124,15 @@ namespace SMESH
 
       myGrid = vtkUnstructuredGrid::New();
 
+      myGeometryFilter = vtkSmartPointer<vtkGeometryFilter>::New();
+      myGeometryFilter->SetInputData(myGrid);
+
+      myTriangleFilter = vtkSmartPointer<vtkTriangleFilter>::New();
+      myTriangleFilter->SetInputConnection(myGeometryFilter->GetOutputPort());
+
       // Create and display actor
       myMapper = vtkDataSetMapper::New();
-      myMapper->SetInputData(myGrid);
+      myMapper->SetInputConnection(myTriangleFilter->GetOutputPort());
 
       myPreviewActor = SALOME_Actor::New();
       myPreviewActor->PickableOff();
@@ -148,7 +158,7 @@ namespace SMESH
 
       // Orientation of faces
       myFaceOrientationFilter = SMESH_FaceOrientationFilter::New();
-      myFaceOrientationFilter->SetInputData(myGrid);
+      myFaceOrientationFilter->SetInputConnection(myTriangleFilter->GetOutputPort());
 
       myFaceOrientationDataMapper = vtkPolyDataMapper::New();
       myFaceOrientationDataMapper->SetInputConnection(myFaceOrientationFilter->GetOutputPort());
@@ -209,6 +219,9 @@ namespace SMESH
       myGrid->InsertNextCell(theType,anIds);
       anIds->Delete();
 
+      myGeometryFilter->Update();
+      myTriangleFilter->Update();
+
       myGrid->Modified();
 
       SetVisibility(true, theActor->GetFacesOriented(), false);

src/SMESHGUI/SMESHGUI_Meshio.cxx

@@ -142,12 +142,13 @@ const QStringList& SMESHGUI_Meshio::GetExportFileFilter()
 #endif
     "FLAC3D (*.f3grid)",
     "Gmsh 2.2 (*.msh)",
-    "Gmsh 4.0, and 4.1 (*.msh)",
+    "Gmsh 4.0 (*.msh)",
+    //"Gmsh 4.1 (*.msh)",
     "H5M (*.h5m)",
     "Kratos/MDPA (*.mdpa)",
     "MED/Salome (*.med)",
     "Medit (*.mesh *.meshb)",
-    "Nastran (*.bdf *fem *.nas)",
+    "Nastran (*.bdf *.fem *.nas)",
     "Netgen(*.vol *.vol.gz)",
     "OBJ (*.obj)",
     "OFF (*.off)",

src/SMESHGUI/SMESHGUI_Operations.h

@@ -94,6 +94,7 @@ namespace SMESHOp {
     OpWhatIs                 = 2101,   // MENU MESH  - MESH ELEMENT INFORMATION
     OpStdInfo                = 2102,   // MENU MESH  - MESH STANDARD INFORMATION
     OpFindElementByPoint     = 2103,   // MENU MESH  - FIND ELEMENT BY POINT
+    OpReloadFromFile         = 2104,   // MENU MESH  - RELOAD MESH FROM FILE
     OpUpdate                 = 2200,   // POPUP MENU - UPDATE
     // Controls -----------------------//--------------------------------
     OpFreeNode               = 3000,   // MENU CONTROLS - FREE NODES

src/SMESHGUI/SMESH_msg_en.ts

@@ -272,6 +272,10 @@
         <source>MEN_ADV_INFO</source>
         <translation>Mesh Information</translation>
     </message>
+    <message>
+        <source>MEN_RELOAD_FROM_FILE</source>
+        <translation>Reload from file</translation>
+    </message>
     <message>
         <source>MEN_ALL</source>
         <translation>All</translation>
@@ -3332,6 +3336,10 @@ Use Display Entity menu command to show them.
         <source>STB_ADV_INFO</source>
         <translation>Show base information about the mesh object</translation>
     </message>
+    <message>
+        <source>STB_RELOAD_FROM_FILE</source>
+        <translation>Reload original mesh from file</translation>
+    </message>
     <message>
         <source>STB_ALL</source>
         <translation>All</translation>
@@ -4072,6 +4080,10 @@ Use Display Entity menu command to show them.
         <source>TOP_ADV_INFO</source>
         <translation>Mesh Information</translation>
     </message>
+    <message>
+        <source>TOP_RELOAD_FROM_FILE</source>
+        <translation>Reload from file</translation>
+    </message>
     <message>
         <source>TOP_ALL</source>
         <translation>All</translation>

src/SMESHGUI/SMESH_msg_fr.ts

@@ -268,6 +268,10 @@
         <source>MEN_ADV_INFO</source>
         <translation>Informations sur le maillage</translation>
     </message>
+    <message>
+        <source>MEN_RELOAD_FROM_FILE</source>
+        <translation>Reload from file</translation>
+    </message>
     <message>
         <source>MEN_ALL</source>
         <translation>Tous</translation>
@@ -3331,6 +3335,10 @@ Utilisez le menu &quot;Visualiser une entité&quot; pour les afficher.
         <source>STB_ALL</source>
         <translation>Tous</translation>
     </message>
+    <message>
+        <source>STB_RELOAD_FROM_FILE</source>
+        <translation>Reload original mesh from file</translation>
+    </message>
     <message>
         <source>STB_AREA</source>
         <translation>Aire</translation>
@@ -4067,6 +4075,10 @@ Utilisez le menu &quot;Visualiser une entité&quot; pour les afficher.
         <source>TOP_ADV_INFO</source>
         <translation>Informations sur le maillage</translation>
     </message>
+    <message>
+        <source>TOP_RELOAD_FROM_FILE</source>
+        <translation>Reload from file</translation>
+    </message>
     <message>
         <source>TOP_ALL</source>
         <translation>Tous</translation>

src/SMESHGUI/SMESH_msg_ja.ts

@@ -251,6 +251,10 @@
       <source>MEN_ADV_INFO</source>
       <translation>メッシュに関する情報</translation>
     </message>
+    <message>
+        <source>MEN_RELOAD_FROM_FILE</source>
+        <translation>Reload from file</translation>
+    </message>
     <message>
       <source>MEN_ALL</source>
       <translation>すべて</translation>
@@ -2979,6 +2983,10 @@ pip install meshio[all]</translation>
       <source>STB_ADV_INFO</source>
       <translation>メッシュ上の基本的な情報を得る</translation>
     </message>
+    <message>
+        <source>STB_RELOAD_FROM_FILE</source>
+        <translation>Reload original mesh from file</translation>
+    </message>
     <message>
       <source>STB_ALL</source>
       <translation>すべて</translation>
@@ -3663,6 +3671,10 @@ pip install meshio[all]</translation>
       <source>TOP_ADV_INFO</source>
       <translation>メッシュに関する情報</translation>
     </message>
+    <message>
+        <source>TOP_RELOAD_FROM_FILE</source>
+        <translation>Reload from file</translation>
+    </message>
     <message>
       <source>TOP_ALL</source>
       <translation>すべて</translation>

src/SMESHUtils/CMakeLists.txt

@@ -23,6 +23,7 @@ INCLUDE_DIRECTORIES(
   ${KERNEL_INCLUDE_DIRS}
   ${OpenCASCADE_INCLUDE_DIR}
   ${Boost_INCLUDE_DIRS}
+  ${SALOMEBOOTSTRAP_INCLUDE_DIRS}
   ${PROJECT_SOURCE_DIR}/src/SMDS
 )
 
@@ -38,6 +39,7 @@ SET(_link_LIBRARIES
    ${OpenCASCADE_FoundationClasses_LIBRARIES}
    ${OpenCASCADE_ModelingData_LIBRARIES}
    ${Boost_LIBRARIES}
+   SALOMEException
    SMDS
 )
 

src/SMESHUtils/SMESH_RegularGrid.cxx

@@ -335,6 +335,71 @@ std::tuple<int,int,int,int> SMESH_RegularGrid::getFaceLimits( const FaceType fac
   }
 }
 
+// Find the left bottom corner of the box face from the 4 points and return its index
+// Auxiliary for the GetCommontInterface function
+static int findLeftBottomCorner(const std::shared_ptr<gp_Pnt>& V0, 
+                                const std::shared_ptr<gp_Pnt>& V1,
+                                const std::shared_ptr<gp_Pnt>& V2,
+                                const std::shared_ptr<gp_Pnt>& V3)
+{
+  // Find the left bottom corner of the box face from the 4 points
+  auto isLeftBottom = [](const std::shared_ptr<gp_Pnt>& a, const std::shared_ptr<gp_Pnt>& b) {
+    return (a->X() < b->X()) || 
+           (a->X() == b->X() && a->Y() < b->Y()) || 
+           (a->X() == b->X() && a->Y() == b->Y() && a->Z() < b->Z());
+  };
+
+  int index = 0;
+
+  if (isLeftBottom(V1, V0)) {
+    index = 1;
+  }
+  if (isLeftBottom(V2, (index == 1 ? V1 : V0))) {
+    index = 2;
+  }
+  if (isLeftBottom(V3, (index == 2 ? V2 : (index == 1 ? V1 : V0)))) {
+    index = 3;
+  }
+
+  return index;
+}
+
+// This function transforms the index limits of a face based on the position of the left bottom corner.
+// It ensures that the face limits are correctly oriented by swapping indices if necessary.
+// Auxiliary for the GetCommontInterface function
+void SMESH_RegularGrid::transformIndexLimits(FaceType face,
+                                             std::vector<int>& faceLimits,
+                                             const std::shared_ptr<gp_Pnt>& V0, 
+                                             const std::shared_ptr<gp_Pnt>& V1,
+                                             const std::shared_ptr<gp_Pnt>& V2, 
+                                             const std::shared_ptr<gp_Pnt>& V3)
+{
+  int indexOfLeftBottom = findLeftBottomCorner(V0, V1, V2, V3);
+  if(indexOfLeftBottom == 0)
+    return;
+
+  if(indexOfLeftBottom == 1)
+  {
+    std::swap(faceLimits[0],faceLimits[3]);
+    std::swap(faceLimits[1],faceLimits[4]);
+    std::swap(faceLimits[2],faceLimits[5]);
+    return;
+  }
+
+  if(face == FaceType::B_BOTTOM || face == FaceType::B_TOP)
+  {
+    indexOfLeftBottom == 2? std::swap(faceLimits[0],faceLimits[3]) : std::swap(faceLimits[1],faceLimits[4]);
+  }
+  else if(face == FaceType::B_RIGHT || face == FaceType::B_LEFT)
+  {
+    indexOfLeftBottom == 2? std::swap(faceLimits[1],faceLimits[4]) : std::swap(faceLimits[2],faceLimits[5]);
+  }
+  else if(face == FaceType::B_BACK || face == FaceType::B_FRONT)
+  {
+    indexOfLeftBottom == 2? std::swap(faceLimits[0],faceLimits[3]) : std::swap(faceLimits[2],faceLimits[5]);
+  }
+}
+
 void SMESH_RegularGrid::GetCommontInterface( FaceType face, SMESH_RegularGrid * grid, std::vector<int>& interface )
 {
   const double tol = Precision::Confusion(); /*confusion is 1e-7, the recommended tolerance to find coincident points in 3D*/
@@ -368,12 +433,11 @@ void SMESH_RegularGrid::GetCommontInterface( FaceType face, SMESH_RegularGrid *
     auto neighboorVertex = grid->getFaceLimits( gridFace );
 
     auto v4 = grid->GetNode(std::get<0>(neighboorVertex));
-    auto v5 = grid->GetNode(std::get<1>(neighboorVertex));  
+    auto v5 = grid->GetNode(std::get<1>(neighboorVertex));
     auto v6 = grid->GetNode(std::get<2>(neighboorVertex));
     auto v7 = grid->GetNode(std::get<3>(neighboorVertex));
 
-
-    std::vector<bool> trueTable({ v0->IsEqual(*v4,tol),v0->IsEqual(*v5,tol),v0->IsEqual(*v6,tol),v0->IsEqual(*v7,tol),                                  
+    std::vector<bool> trueTable({ v0->IsEqual(*v4,tol),v0->IsEqual(*v5,tol),v0->IsEqual(*v6,tol),v0->IsEqual(*v7,tol),
                                   v1->IsEqual(*v4,tol),v1->IsEqual(*v5,tol),v1->IsEqual(*v6,tol),v1->IsEqual(*v7,tol),
                                   v2->IsEqual(*v4,tol),v2->IsEqual(*v5,tol),v2->IsEqual(*v6,tol),v2->IsEqual(*v7,tol),
                                   v3->IsEqual(*v4,tol),v3->IsEqual(*v5,tol),v3->IsEqual(*v6,tol),v3->IsEqual(*v7,tol)});
@@ -389,7 +453,7 @@ void SMESH_RegularGrid::GetCommontInterface( FaceType face, SMESH_RegularGrid *
          trueCorner[3] )   /*Face to Face interface 100% conform*/
     {
       interfaceRange = this->getFaceIndexLimits<int>( face );
-      interfaceDonor = grid->getFaceIndexLimits<int>( gridFace );      
+      interfaceDonor = grid->getFaceIndexLimits<int>( gridFace );
     }
     else if ( trueCorner[0] || trueCorner[1]  ||
               trueCorner[2] || trueCorner[3] ) /*Partial Face to Face. Only one intersection then all the other 3 vertex are: 2 in the edges 1 inside the face*/
@@ -401,7 +465,7 @@ void SMESH_RegularGrid::GetCommontInterface( FaceType face, SMESH_RegularGrid *
         auto nodeToSearch = trueCorner[0] ? v1 : trueCorner[2] ? v3 : trueCorner[1] ? v0 : v2;
 
         grid->foreachNodeOnFace( gridFace, [&] (const std::shared_ptr<gp_Pnt> sidePoint, const int nodeIndex) 
-        {       
+        {
           if ( nodeToSearch->IsEqual( *sidePoint, tol ) ) 
           {
             interfaceRange = this->getFaceIndexLimits<int>( face );
@@ -415,7 +479,7 @@ void SMESH_RegularGrid::GetCommontInterface( FaceType face, SMESH_RegularGrid *
                               grid->getFaceIndexLimits( startIndex, nodeIndex ) : 
                                 grid->getFaceIndexLimits( nodeIndex, startIndex );
             return;
-          }          
+          }
         });
       }
       // 2) face > gridEdge
@@ -437,23 +501,27 @@ void SMESH_RegularGrid::GetCommontInterface( FaceType face, SMESH_RegularGrid *
                                 this->getFaceIndexLimits( startIndex, nodeIndex) : 
                                   this->getFaceIndexLimits( nodeIndex, startIndex );
 
-            interfaceDonor = grid->getFaceIndexLimits<int>( gridFace );  
+            interfaceDonor = grid->getFaceIndexLimits<int>( gridFace );
             return;
-          }          
+          }
         });
-      }      
+      }
     }
 
     if ( !interfaceRange.empty() && !interfaceDonor.empty() )
-    {      
+    {
+      // Transform the index limits of the face based on the position of the left bottom corner
+      transformIndexLimits(face, interfaceRange, v0, v1, v2, v3);
+      transformIndexLimits(gridFace, interfaceDonor, v4, v5, v6, v7);
+
       tranformationRange  = this->computeTransformation( face, gridFace, interfaceRange, interfaceDonor );
       tranformationDonor  = grid->computeTransformation( gridFace, face, interfaceDonor, interfaceRange );
       interface           = std::vector<int>{(int)face}+interfaceRange+interfaceDonor+tranformationRange;
       auto dualInterface  = std::vector<int>{(int)face}+interfaceDonor+interfaceRange+tranformationDonor;
       this->setInterface( face, grid->id(), interface );
       grid->setInterface( gridFace, this->id(), dualInterface );
-    }    
-  });    
+    }
+  });
 }
 
 void SMESH_RegularGrid::GetCommontInterface( EdgeType edge, SMESH_RegularGrid * grid, std::vector<int>& interface )

src/SMESHUtils/SMESH_RegularGrid.hxx

@@ -353,6 +353,16 @@ namespace SMESHUtils
     // \return the index of all the nodes in the face 
     std::vector<int> nodesOfFace( SMESH_RegularGrid::FaceType face ) const;
 
+    // This function transforms the index limits of a face based on the position of the left bottom corner.
+    // It ensures that the face limits are correctly oriented by swapping indices if necessary.
+    // Auxiliary for the GetCommontInterface function
+    void transformIndexLimits(FaceType face, 
+                              std::vector<int>& faceLimits, 
+                              const std::shared_ptr<gp_Pnt>& V0, 
+                              const std::shared_ptr<gp_Pnt>& V1,
+                              const std::shared_ptr<gp_Pnt>& V2, 
+                              const std::shared_ptr<gp_Pnt>& V3);
+
   private:
     int myId;
     int mnx,mny,mnz,mns;

src/SMESHUtils/SMESH_RegularGridTemplate.hxx

@@ -80,6 +80,40 @@ namespace SMESHUtils
     }
   }
 
+  /*    Define vertex numeration convention for vertex and face index limits!
+  *
+  *               V7                           V6
+  *                +--------+--------+--------+
+  *               /        /        /        /|
+  *              +--------+--------+--------+ |
+  *             /        /        /        /| |
+  *            +--------+--------+--------+ | +
+  *           /        /        /     V5 /| |/|
+  *       V4 +--------+--------+--------+ | + |
+  *          |     |                    | |/| |
+  *          |     |                    | + | +
+  *          |     |                    |/| |/|
+  *          +     |                    | + |
+  *          |     |                    | |/| |
+  *          |     +-------+----------+-| + | + V2
+  *          |    / V3                  |/| |/
+  *          +   /                      | +
+  *          |  /                       | |/
+  *          | /                        | +
+  *          |/                         |/
+  *          +--------+--------+--------+
+  *         V0                          V1
+  * 
+  * Canonical cartesian axis orientation
+  * 
+  *         ^   ^ j (or y)
+  * k (or z)|  /
+  *         | /
+  *         |/
+  *         +------> i (or x)
+  * 
+  * 
+  */
   template<typename T>
   std::vector<T> SMESH_RegularGrid::getFaceIndexLimits( const FaceType face ) const
   {

src/SMESH_I/SMESH_2smeshpy.cxx

@@ -1013,6 +1013,7 @@ void _pyGen::Process( const Handle(_pyCommand)& theCommand )
 {
   // there are methods to convert:
   // CreateMesh( shape )
+  // ReloadMesh( shape )
   // Concatenate( [mesh1, ...], ... )
   // CreateHypothesis( theHypType, theLibName )
   // Compute( mesh, geom )
@@ -1067,6 +1068,13 @@ void _pyGen::Process( const Handle(_pyCommand)& theCommand )
     Handle(_pyMesh) mesh = new _pyMesh( theCommand, entries.front() );
     AddObject( mesh );
   }
+  if (method == "ReloadMeshFromFile")
+  {
+    std::cout << "WhatIS" << std::endl;
+    Handle(_pyMesh) mesh = new _pyMesh(theCommand, theCommand->GetResultValue());
+    AddObject(mesh);
+    return;
+  }
 
   // CreateHypothesis()
   if ( method == "CreateHypothesis" )

src/SMESH_I/SMESH_Gen_i.cxx

@@ -153,6 +153,7 @@
 #include <cstdio>
 #include <cstdlib>
 #include <memory>
+#include <QStringList>
 
 #include <boost/archive/text_oarchive.hpp>
 #include <boost/serialization/list.hpp>
@@ -1330,6 +1331,61 @@ SMESH::SMESH_Mesh_ptr SMESH_Gen_i::CreateEmptyMesh()
   return mesh._retn();
 }
 
+SMESH::SMESH_Mesh_ptr SMESH_Gen_i::ReloadMeshFromFile(SMESH::SMESH_Mesh_ptr theMesh)
+{
+  SMESH::mesh_array_var aMeshes = new SMESH::mesh_array;
+
+  SMESH::MedFileInfo* aMedInfo = theMesh->GetMEDFileInfo();
+
+  // Get file path and re-import mesh
+  QString aPath = QString(aMedInfo->fileName);
+  QStringList aSplit = aPath.split('.');
+  QStringList anEntryList;
+
+  auto aStudy = getStudyServant();
+
+  SMESH::SMESH_Mesh_ptr aNewMesh;
+
+  SMESH_Mesh* aMesh = reinterpret_cast<SMESH_Mesh*> (theMesh->GetMeshPtr());
+  aMesh->GetMeshDS()->ClearMesh();
+
+  if (aSplit.last() == "cgns")
+  {
+    SMESH::DriverMED_ReadStatus res;
+    aMeshes = ReloadMeshesFromCGNS(aPath.toUtf8().constData(), theMesh, res);
+
+    aNewMesh = aMeshes[0];
+
+  }
+  else if (aSplit.last().contains("stl", Qt::CaseSensitivity::CaseInsensitive))
+  {
+    aNewMesh = ReloadMeshesFromSTL(aPath.toUtf8().constData(), theMesh);
+  }
+  else if (aSplit.last().contains("unv", Qt::CaseSensitivity::CaseInsensitive))
+  {
+    aNewMesh = ReloadMeshesFromUNV(aPath.toUtf8().constData(), theMesh);
+  }
+  else if (aSplit.last().contains("mesh", Qt::CaseSensitivity::CaseInsensitive))
+  {
+    SMESH::ComputeError_var res;
+    aNewMesh = ReloadMeshesFromGMF(aPath.toUtf8().constData(), theMesh, true, res.out());
+  }
+  else if (aSplit.last().contains("med", Qt::CaseSensitivity::CaseInsensitive))
+  {
+    SMESH::DriverMED_ReadStatus res;
+    aMeshes = ReloadMeshesFromMED(aPath.toUtf8().constData(), theMesh, res);
+
+    aNewMesh = aMeshes[0];
+  }
+  else
+  {
+    // MeshIO
+  }
+
+  theMesh = SMESH::SMESH_Mesh::_duplicate(aNewMesh);
+  return theMesh;
+}
+
 namespace
 {
   //================================================================================
@@ -1401,12 +1457,15 @@ SMESH::SMESH_Mesh_ptr SMESH_Gen_i::CreateMeshesFromUNV( const char* theFileName
   checkFileReadable( theFileName );
 
   SMESH::SMESH_Mesh_var aMesh = createMesh();
+
   string aFileName;
+
   // publish mesh in the study
   if ( CanPublishInStudy( aMesh ) ) {
     SALOMEDS::StudyBuilder_var aStudyBuilder = getStudyServant()->NewBuilder();
     aStudyBuilder->NewCommand();  // There is a transaction
-    SALOMEDS::SObject_wrap aSO = PublishMesh( aMesh.in(), aFileName.c_str() );
+    SALOMEDS::SObject_wrap aSO =
+      PublishMesh( aMesh.in(), aFileName.c_str(), "ICON_SMESH_TREE_MESH_IMPORTED" );
     aStudyBuilder->CommitCommand();
     if ( !aSO->_is_nil() ) {
       // Update Python script
@@ -1425,6 +1484,36 @@ SMESH::SMESH_Mesh_ptr SMESH_Gen_i::CreateMeshesFromUNV( const char* theFileName
   return aMesh._retn();
 }
 
+SMESH::SMESH_Mesh_ptr SMESH_Gen_i::ReloadMeshesFromUNV(const char* theFileName, SMESH::SMESH_Mesh_ptr sourceMesh)
+{
+  Unexpect aCatch(SALOME_SalomeException);
+
+  checkFileReadable(theFileName);
+
+  string aFileName;
+  // publish mesh in the study
+  if (CanPublishInStudy(sourceMesh)) {
+    SALOMEDS::StudyBuilder_var aStudyBuilder = getStudyServant()->NewBuilder();
+    aStudyBuilder->NewCommand();  // There is a transaction
+    SALOMEDS::SObject_wrap aSO = PublishMesh(sourceMesh, aFileName.c_str());
+    aStudyBuilder->CommitCommand();
+    if (!aSO->_is_nil()) {
+      // Update Python script
+      TPythonDump(this) << aSO << " = " << this << ".ReloadMeshesFromUNV(r'" << theFileName << "')";
+    }
+  }
+
+  SMESH_Mesh_i* aServant = dynamic_cast<SMESH_Mesh_i*>(GetServant(sourceMesh).in());
+  ASSERT(aServant);
+  aServant->ImportUNVFile(theFileName);
+
+  // Dump creation of groups
+  SMESH::ListOfGroups_var groups = aServant->GetGroups();
+
+  aServant->GetImpl().GetMeshDS()->Modified();
+  return sourceMesh;
+}
+
 //=============================================================================
 /*!
  *  SMESH_Gen_i::CreateMeshFromMED
@@ -1438,13 +1527,6 @@ SMESH::mesh_array* SMESH_Gen_i::CreateMeshesFromMED( const char*
 {
   checkFileReadable( theFileName );
 
-#ifdef WIN32
-  char bname[ _MAX_FNAME ];
-  _splitpath( theFileName, NULL, NULL, bname, NULL );
-  string aFileName = bname;
-#else
-  string aFileName = basename( const_cast<char *>( theFileName ));
-#endif
   // Retrieve mesh names from the file
   DriverMED_R_SMESHDS_Mesh myReader;
   myReader.SetFile( theFileName );
@@ -1481,7 +1563,7 @@ SMESH::mesh_array* SMESH_Gen_i::CreateMeshesFromMED( const char*
         // publish mesh in the study
         SALOMEDS::SObject_wrap aSO;
         if ( CanPublishInStudy( mesh ) )
-          aSO = PublishMesh( mesh.in(), meshName.c_str() );
+          aSO = PublishMesh( mesh.in(), meshName.c_str(), "ICON_SMESH_TREE_MESH_IMPORTED" );
 
         // Python Dump
         if ( !aSO->_is_nil() ) {
@@ -1515,6 +1597,92 @@ SMESH::mesh_array* SMESH_Gen_i::CreateMeshesFromMED( const char*
   return aResult._retn();
 }
 
+SMESH::mesh_array* SMESH_Gen_i::ReloadMeshesFromMED(const char* theFileName, SMESH::SMESH_Mesh_ptr sourceMesh, SMESH::DriverMED_ReadStatus& theStatus)
+{
+  SMESH::ListOfGroups anOldGroups = *sourceMesh->GetGroups();
+  checkFileReadable(theFileName);
+
+#ifdef WIN32
+  char bname[_MAX_FNAME];
+  _splitpath(theFileName, NULL, NULL, bname, NULL);
+  string aFileName = bname;
+#else
+  string aFileName = basename(const_cast<char*>(theFileName));
+#endif
+  // Retrieve mesh names from the file
+  DriverMED_R_SMESHDS_Mesh myReader;
+  myReader.SetFile(theFileName);
+  myReader.SetMeshId(-1);
+  Driver_Mesh::Status aStatus;
+  list<string> aNames = myReader.GetMeshNames(aStatus);
+  SMESH::mesh_array_var aResult = new SMESH::mesh_array();
+  theStatus = (SMESH::DriverMED_ReadStatus)aStatus;
+
+  { // open a new scope to make aPythonDump die before PythonDump in SMESH_Mesh::GetGroups()
+
+    // Python Dump
+    TPythonDump aPythonDump(this);
+    aPythonDump << "([";
+
+    if (theStatus == SMESH::DRS_OK)
+    {
+      SALOMEDS::StudyBuilder_var aStudyBuilder;
+      aStudyBuilder = getStudyServant()->NewBuilder();
+      aStudyBuilder->NewCommand();  // There is a transaction
+
+      aResult->length(aNames.size());
+      int i = 0;
+
+      // Iterate through all meshes and create mesh objects
+      for (const std::string& meshName : aNames)
+      {
+        // Python Dump
+        if (i > 0) aPythonDump << ", ";
+
+        // publish mesh in the study
+        SALOMEDS::SObject_wrap aSO;
+        if (CanPublishInStudy(sourceMesh))
+          aSO = PublishMesh(sourceMesh, meshName.c_str());
+
+        // Python Dump
+        if (!aSO->_is_nil()) {
+          aPythonDump << aSO;
+        }
+        else {
+          aPythonDump << "mesh_" << i;
+        }
+
+        // Read mesh data (groups are published automatically by ImportMEDFile())
+        SMESH_Mesh_i* meshServant = dynamic_cast<SMESH_Mesh_i*>(GetServant(sourceMesh).in());
+        ASSERT(meshServant);
+        SMESH::DriverMED_ReadStatus status1 =
+          meshServant->ImportMEDFile(theFileName, meshName.c_str());
+        if (status1 > theStatus)
+          theStatus = status1;
+
+        for (unsigned int i = 0; i < anOldGroups.length(); ++i)
+        {
+          //auto X = anOldGroups[i];
+          sourceMesh->RemoveGroup(anOldGroups[i]);
+        }
+        aResult[i++] = SMESH::SMESH_Mesh::_duplicate(sourceMesh);
+
+        meshServant->GetImpl().GetMeshDS()->Modified();
+      }
+      if (!aStudyBuilder->_is_nil())
+        aStudyBuilder->CommitCommand();
+    }
+
+    // Update Python script
+    aPythonDump << "], status) = " << this << ".ReloadMeshesFromMED( r'" << theFileName << "' )";
+  }
+  // Dump creation of groups
+  for (CORBA::ULong i = 0; i < aResult->length(); ++i)
+    SMESH::ListOfGroups_var groups = aResult[i]->GetGroups();
+
+  return aResult._retn();
+}
+
 //=============================================================================
 /*!
  *  SMESH_Gen_i::CreateMeshFromSTL
@@ -1529,7 +1697,7 @@ SMESH::SMESH_Mesh_ptr SMESH_Gen_i::CreateMeshesFromSTL( const char* theFileName
   checkFileReadable( theFileName );
 
   SMESH::SMESH_Mesh_var aMesh = createMesh();
-  //string aFileName;
+
 #ifdef WIN32
   char bname[ _MAX_FNAME ];
   _splitpath( theFileName, NULL, NULL, bname, NULL );
@@ -1537,11 +1705,13 @@ SMESH::SMESH_Mesh_ptr SMESH_Gen_i::CreateMeshesFromSTL( const char* theFileName
 #else
   string aFileName = basename( const_cast<char *>(theFileName) );
 #endif
+
   // publish mesh in the study
   if ( CanPublishInStudy( aMesh ) ) {
     SALOMEDS::StudyBuilder_var aStudyBuilder = getStudyServant()->NewBuilder();
     aStudyBuilder->NewCommand();  // There is a transaction
-    SALOMEDS::SObject_wrap aSO = PublishInStudy( SALOMEDS::SObject::_nil(), aMesh.in(), aFileName.c_str() );
+    SALOMEDS::SObject_wrap aSO =
+      PublishMesh( aMesh.in(), aFileName.c_str(), "ICON_SMESH_TREE_MESH_IMPORTED" );
     aStudyBuilder->CommitCommand();
     if ( !aSO->_is_nil() ) {
       // Update Python script
@@ -1556,6 +1726,37 @@ SMESH::SMESH_Mesh_ptr SMESH_Gen_i::CreateMeshesFromSTL( const char* theFileName
   return aMesh._retn();
 }
 
+SMESH::SMESH_Mesh_ptr SMESH_Gen_i::ReloadMeshesFromSTL(const char* theFileName, SMESH::SMESH_Mesh_ptr sourceMesh)
+{
+  Unexpect aCatch(SALOME_SalomeException);
+  checkFileReadable(theFileName);
+
+#ifdef WIN32
+  char bname[_MAX_FNAME];
+  _splitpath(theFileName, NULL, NULL, bname, NULL);
+  string aFileName = bname;
+#else
+  string aFileName = basename(const_cast<char*>(theFileName));
+#endif
+  // publish mesh in the study
+  if (CanPublishInStudy(sourceMesh)) {
+    SALOMEDS::StudyBuilder_var aStudyBuilder = getStudyServant()->NewBuilder();
+    aStudyBuilder->NewCommand();  // There is a transaction
+    SALOMEDS::SObject_wrap aSO = PublishInStudy(SALOMEDS::SObject::_nil(), sourceMesh, aFileName.c_str());
+    aStudyBuilder->CommitCommand();
+    if (!aSO->_is_nil()) {
+      // Update Python script
+      TPythonDump(this) << aSO << " = " << this << ".ReloadMeshFromSTL(r'" << theFileName << "')";
+    }
+  }
+
+  SMESH_Mesh_i* aServant = dynamic_cast<SMESH_Mesh_i*>(GetServant(sourceMesh).in());
+  ASSERT(aServant);
+  aServant->ImportSTLFile(theFileName);
+  aServant->GetImpl().GetMeshDS()->Modified();
+  return sourceMesh;
+}
+
 //================================================================================
 /*!
  * \brief Create meshes and import data from the CGSN file
@@ -1616,7 +1817,7 @@ SMESH::mesh_array* SMESH_Gen_i::CreateMeshesFromCGNS( const char*
         // publish mesh in the study
         SALOMEDS::SObject_wrap aSO;
         if ( CanPublishInStudy( mesh ) )
-          aSO = PublishMesh( mesh.in(), meshName.c_str() );
+          aSO = PublishMesh( mesh.in(), meshName.c_str(), "ICON_SMESH_TREE_MESH_IMPORTED" );
 
         // Python Dump
         if ( !aSO->_is_nil() ) {
@@ -1641,6 +1842,85 @@ SMESH::mesh_array* SMESH_Gen_i::CreateMeshesFromCGNS( const char*
   return aResult._retn();
 }
 
+SMESH::mesh_array* SMESH_Gen_i::ReloadMeshesFromCGNS(const char*                  theFileName,
+                                                     SMESH::SMESH_Mesh_ptr        sourceMesh,
+                                                     SMESH::DriverMED_ReadStatus& theStatus)
+{
+  Unexpect aCatch(SALOME_SalomeException);
+  checkFileReadable(theFileName);
+
+  SMESH::mesh_array_var aResult = new SMESH::mesh_array();
+
+#ifdef WITH_CGNS
+  // Retrieve nb meshes from the file
+  DriverCGNS_Read myReader;
+  myReader.SetFile(theFileName);
+  Driver_Mesh::Status aStatus;
+  int nbMeshes = myReader.GetNbMeshes(aStatus);
+  theStatus = (SMESH::DriverMED_ReadStatus)aStatus;
+
+  aResult->length(nbMeshes);
+
+  { // open a new scope to make aPythonDump die before PythonDump in SMESH_Mesh::GetGroups()
+
+    // Python Dump
+    TPythonDump aPythonDump(this);
+    aPythonDump << "([";
+
+    if (theStatus == SMESH::DRS_OK)
+    {
+      SALOMEDS::StudyBuilder_var aStudyBuilder = getStudyServant()->NewBuilder();
+      aStudyBuilder->NewCommand();  // There is a transaction
+
+      int i = 0;
+
+      // Iterate through all meshes and create mesh objects
+      for (; i < nbMeshes; ++i)
+      {
+        // Python Dump
+        if (i > 0) aPythonDump << ", ";
+
+        // create mesh
+        aResult[i] = SMESH::SMESH_Mesh::_duplicate(sourceMesh);
+
+        // Read mesh data (groups are published automatically by ImportMEDFile())
+        SMESH_Mesh_i* meshServant = dynamic_cast<SMESH_Mesh_i*>(GetServant(sourceMesh).in());
+        ASSERT(meshServant);
+        string meshName;
+        SMESH::DriverMED_ReadStatus status1 =
+          meshServant->ImportCGNSFile(theFileName, i, meshName);
+        if (status1 > theStatus)
+          theStatus = status1;
+
+        meshServant->GetImpl().GetMeshDS()->Modified();
+        // publish mesh in the study
+        SALOMEDS::SObject_wrap aSO;
+        if (CanPublishInStudy(sourceMesh))
+          aSO = PublishMesh(sourceMesh, meshName.c_str());
+
+        // Python Dump
+        if (!aSO->_is_nil()) {
+          aPythonDump << aSO;
+        }
+        else {
+          aPythonDump << "mesh_" << i;
+        }
+      }
+      aStudyBuilder->CommitCommand();
+    }
+
+    aPythonDump << "], status) = " << this << ".ReloadMeshesFromCGNS(r'" << theFileName << "')";
+  }
+  // Dump creation of groups
+  for (CORBA::ULong i = 0; i < aResult->length(); ++i)
+    SMESH::ListOfGroups_var groups = aResult[i]->GetGroups();
+#else
+  THROW_SALOME_CORBA_EXCEPTION("CGNS library is unavailable", SALOME::INTERNAL_ERROR);
+#endif
+
+  return aResult._retn();
+}
+
 //================================================================================
 /*!
  * \brief Create a mesh and import data from a GMF file
@@ -1667,7 +1947,8 @@ SMESH_Gen_i::CreateMeshesFromGMF( const char*             theFileName,
   if ( CanPublishInStudy( aMesh ) ) {
     SALOMEDS::StudyBuilder_var aStudyBuilder = getStudyServant()->NewBuilder();
     aStudyBuilder->NewCommand();  // There is a transaction
-    SALOMEDS::SObject_wrap aSO = PublishInStudy( SALOMEDS::SObject::_nil(), aMesh.in(), aFileName.c_str() );
+    SALOMEDS::SObject_wrap aSO =
+      PublishMesh( aMesh.in(), aFileName.c_str(), "ICON_SMESH_TREE_MESH_IMPORTED" );
     aStudyBuilder->CommitCommand();
     if ( !aSO->_is_nil() ) {
       // Update Python script
@@ -1683,6 +1964,38 @@ SMESH_Gen_i::CreateMeshesFromGMF( const char*             theFileName,
   return aMesh._retn();
 }
 
+SMESH::SMESH_Mesh_ptr SMESH_Gen_i::ReloadMeshesFromGMF(const char* theFileName, SMESH::SMESH_Mesh_ptr sourceMesh, CORBA::Boolean theMakeRequiredGroups, SMESH::ComputeError_out theError)
+{
+  Unexpect aCatch(SALOME_SalomeException);
+  checkFileReadable(theFileName);
+
+#ifdef WIN32
+  char bname[_MAX_FNAME];
+  _splitpath(theFileName, NULL, NULL, bname, NULL);
+  string aFileName = bname;
+#else
+  string aFileName = basename(const_cast<char*>(theFileName));
+#endif
+  // publish mesh in the study
+  if (CanPublishInStudy(sourceMesh)) {
+    SALOMEDS::StudyBuilder_var aStudyBuilder = getStudyServant()->NewBuilder();
+    aStudyBuilder->NewCommand();  // There is a transaction
+    SALOMEDS::SObject_wrap aSO = PublishInStudy(SALOMEDS::SObject::_nil(), sourceMesh, aFileName.c_str());
+    aStudyBuilder->CommitCommand();
+    if (!aSO->_is_nil()) {
+      // Update Python script
+      TPythonDump(this) << "(" << aSO << ", error) = " << this << ".ReloadMeshesFromGMF(r'"
+        << theFileName << "', "
+        << theMakeRequiredGroups << " )";
+    }
+  }
+  SMESH_Mesh_i* aServant = dynamic_cast<SMESH_Mesh_i*>(GetServant(sourceMesh).in());
+  ASSERT(aServant);
+  theError = aServant->ImportGMFFile(theFileName, theMakeRequiredGroups);
+  aServant->GetImpl().GetMeshDS()->Modified();
+  return sourceMesh;
+}
+
 //================================================================================
 /*!
  * \brief Create a mesh and import data from any file supported by meshio library
@@ -1737,7 +2050,7 @@ SMESH::mesh_array* SMESH_Gen_i::CreateMeshesFromMESHIO(const char* theFileName,
       // publish mesh in the study
       SALOMEDS::SObject_wrap aSO;
       if (CanPublishInStudy(mesh))
-        aSO = PublishMesh(mesh.in(), meshName.c_str());
+        aSO = PublishMesh(mesh.in(), meshName.c_str(), "ICON_SMESH_TREE_MESH_IMPORTED");
 
       // Save SO to use in a python dump
       sobjects.emplace_back(aSO);
@@ -1769,6 +2082,82 @@ SMESH::mesh_array* SMESH_Gen_i::CreateMeshesFromMESHIO(const char* theFileName,
   return aResult._retn();
 }
 
+SMESH::mesh_array* SMESH_Gen_i::ReloadMeshesFromMESHIO(const char* theFileName, SMESH::SMESH_Mesh_ptr sourceMesh, SMESH::DriverMED_ReadStatus& theStatus)
+{
+  Unexpect aCatch(SALOME_SalomeException);
+  checkFileReadable(theFileName);
+
+  MESSAGE("Import part with meshio through an intermediate MED file");
+
+  // Create an object that holds a temp file name and
+  // removes the file when goes out of scope.
+  SMESH_Meshio meshio;
+  const QString tempFileName = meshio.CreateTempFileName(theFileName);
+
+  // Convert temp file into a target one with meshio command
+  meshio.Convert(theFileName, tempFileName);
+
+  // We don't need a python dump from SMESH_Gen_i::CreateMeshesFromMED(), so
+  // we can't use this method as is here. The followed code is an edited part of
+  // copy pasted CreateMeshesFromMED().
+
+  // Retrieve mesh names from the file
+  DriverMED_R_SMESHDS_Mesh myReader;
+  myReader.SetFile(tempFileName.toStdString());
+  myReader.SetMeshId(-1);
+  Driver_Mesh::Status aStatus;
+  list<string> aNames = myReader.GetMeshNames(aStatus);
+  SMESH::mesh_array_var aResult = new SMESH::mesh_array();
+  theStatus = (SMESH::DriverMED_ReadStatus)aStatus;
+
+  if (theStatus == SMESH::DRS_OK)
+  {
+    SALOMEDS::StudyBuilder_var aStudyBuilder;
+    aStudyBuilder = getStudyServant()->NewBuilder();
+    aStudyBuilder->NewCommand();  // There is a transaction
+
+    aResult->length(aNames.size());
+    std::vector<SALOMEDS::SObject_wrap> sobjects;
+    int i = 0;
+
+    // Iterate through all meshes and create mesh objects
+    for (const std::string& meshName : aNames)
+    {
+      // publish mesh in the study
+      SALOMEDS::SObject_wrap aSO;
+      if (CanPublishInStudy(sourceMesh))
+        aSO = PublishMesh(sourceMesh, meshName.c_str());
+
+      // Save SO to use in a python dump
+      sobjects.emplace_back(aSO);
+
+      // Read mesh data (groups are published automatically by ImportMEDFile())
+      SMESH_Mesh_i* meshServant = dynamic_cast<SMESH_Mesh_i*>(GetServant(sourceMesh).in());
+      ASSERT(meshServant);
+      SMESH::DriverMED_ReadStatus status1 =
+        meshServant->ImportMEDFile(tempFileName.toUtf8().data(), meshName.c_str());
+      if (status1 > theStatus)
+        theStatus = status1;
+
+      aResult[i++] = SMESH::SMESH_Mesh::_duplicate(sourceMesh);
+      meshServant->GetImpl().GetMeshDS()->Modified();
+    }
+
+    if (!aStudyBuilder->_is_nil())
+      aStudyBuilder->CommitCommand();
+
+    // Python dump
+    const std::string functionName = std::string(".ReloadMeshesFromMESHIO(r'") + theFileName + "')";
+    functionToPythonDump(this, functionName, sobjects);
+  }
+
+  // Dump creation of groups
+  for (CORBA::ULong i = 0; i < aResult->length(); ++i)
+    SMESH::ListOfGroups_var groups = aResult[i]->GetGroups();
+
+  return aResult._retn();
+}
+
 //=============================================================================
 /*!
  *  SMESH_Gen_i::IsReadyToCompute

src/SMESH_I/SMESH_Gen_i.hxx

@@ -237,29 +237,54 @@ public:
   // Create empty mesh
   SMESH::SMESH_Mesh_ptr CreateEmptyMesh();
 
+  SMESH::SMESH_Mesh_ptr ReloadMeshFromFile(SMESH::SMESH_Mesh_ptr theMesh);
+
   //  Create a mesh and import data from an UNV file
   SMESH::SMESH_Mesh_ptr CreateMeshesFromUNV( const char* theFileName );
 
+  SMESH::SMESH_Mesh_ptr ReloadMeshesFromUNV(const char*           theFileName,
+                                            SMESH::SMESH_Mesh_ptr sourceMesh);
+
   //  Create mesh(es) and import data from MED file
   SMESH::mesh_array* CreateMeshesFromMED( const char* theFileName,
                                           SMESH::DriverMED_ReadStatus& theStatus );
 
+  SMESH::mesh_array* ReloadMeshesFromMED(const char*                  theFileName,
+                                         SMESH::SMESH_Mesh_ptr        sourceMesh,
+                                         SMESH::DriverMED_ReadStatus& theStatus);
+
   //  Create a mesh and import data from a STL file
   SMESH::SMESH_Mesh_ptr CreateMeshesFromSTL( const char* theFileName );
 
+  SMESH::SMESH_Mesh_ptr ReloadMeshesFromSTL(const char*           theFileName,
+                                            SMESH::SMESH_Mesh_ptr sourceMesh);
+
   //  Create mesh(es) and import data from CGNS file
   SMESH::mesh_array* CreateMeshesFromCGNS( const char* theFileName,
                                            SMESH::DriverMED_ReadStatus& theStatus );
 
+  SMESH::mesh_array* ReloadMeshesFromCGNS(const char*                  theFileName,
+                                          SMESH::SMESH_Mesh_ptr        sourceMesh, 
+                                          SMESH::DriverMED_ReadStatus& theStatus);
+
   //  Create a mesh and import data from a GMF file
   SMESH::SMESH_Mesh_ptr CreateMeshesFromGMF( const char*             theFileName,
                                              CORBA::Boolean          theMakeRequiredGroups,
                                              SMESH::ComputeError_out theError);
 
+  SMESH::SMESH_Mesh_ptr ReloadMeshesFromGMF(const char* theFileName,
+                                            SMESH::SMESH_Mesh_ptr        sourceMesh,
+                                            CORBA::Boolean          theMakeRequiredGroups,
+                                            SMESH::ComputeError_out theError);
+
   //  Create a mesh and import data from any file supported by meshio library
   SMESH::mesh_array* CreateMeshesFromMESHIO(const char*             theFileName,
                                             SMESH::DriverMED_ReadStatus& theStatus);
 
+  SMESH::mesh_array* ReloadMeshesFromMESHIO(const char* theFileName,
+                                          SMESH::SMESH_Mesh_ptr        sourceMesh,
+                                          SMESH::DriverMED_ReadStatus& theStatus);
+
   // Create dual mesh of a tetrahedron mesh
   SMESH::SMESH_Mesh_ptr CreateDualMesh(SMESH::SMESH_IDSource_ptr meshPart,
                                        const char*               meshName,
@@ -524,7 +549,8 @@ public:
   // publishing methods
   SALOMEDS::SComponent_ptr PublishComponent();
   SALOMEDS::SObject_ptr PublishMesh (SMESH::SMESH_Mesh_ptr theMesh,
-                                     const char*           theName = 0);
+                                     const char*           theName = 0,
+                                     const char*           thePixMap = 0);
   SALOMEDS::SObject_ptr PublishHypothesis (SMESH::SMESH_Hypothesis_ptr theHyp,
                                            const char*                 theName = 0);
   SALOMEDS::SObject_ptr PublishSubMesh (SMESH::SMESH_Mesh_ptr    theMesh,
@@ -536,6 +562,7 @@ public:
                                       GEOM::GEOM_Object_ptr  theShapeObject,
                                       const char*            theName = 0);
   void UpdateIcons(SMESH::SMESH_Mesh_ptr theMesh);
+  void UpdateGroupIcon(SMESH::SMESH_GroupBase_ptr theGroup);
   void HighLightInvalid(CORBA::Object_ptr theObject, bool isInvalid);
   bool IsInvalid(SALOMEDS::SObject_ptr theObject);
   bool AddHypothesisToShape(SMESH::SMESH_Mesh_ptr       theMesh,

src/SMESH_I/SMESH_Gen_i_1.cxx

@@ -607,7 +607,8 @@ SALOMEDS::SComponent_ptr SMESH_Gen_i::PublishComponent()
 //=======================================================================
 
 SALOMEDS::SObject_ptr SMESH_Gen_i::PublishMesh (SMESH::SMESH_Mesh_ptr theMesh,
-                                                const char*           theName)
+                                                const char*           theName,
+                                                const char*           thePixMap)
 {
   if ( !myIsEnablePublish )
     return SALOMEDS::SObject::_nil();
@@ -631,7 +632,8 @@ SALOMEDS::SObject_ptr SMESH_Gen_i::PublishMesh (SMESH::SMESH_Mesh_ptr theMesh,
     else
       aTag++;
 
-    aMeshSO = publish( theMesh, father, aTag, "ICON_SMESH_TREE_MESH_WARN" );
+    aMeshSO = publish( theMesh, father, aTag,
+                       thePixMap ? thePixMap : "ICON_SMESH_TREE_MESH_WARN" );
     if ( aMeshSO->_is_nil() )
       return aMeshSO._retn();
   }
@@ -977,18 +979,13 @@ void SMESH_Gen_i::UpdateIcons( SMESH::SMESH_Mesh_ptr theMesh )
       if ( idSrc->_is_nil() )
         continue;
 
-      SMESH::SMESH_GroupBase_var     grp = SMESH::SMESH_GroupBase::_narrow( obj );
-      SMESH::SMESH_GroupOnFilter_var gof = SMESH::SMESH_GroupOnFilter::_narrow( obj );
-      const bool         isGroup = !grp->_is_nil();
-      const bool isGroupOnFilter = !gof->_is_nil();
-
       bool isEmpty = ( mesh_i->NbNodes() == 0 );
       if ( !isEmpty )
       {
-        if ( isGroupOnFilter ) // GetTypes() can be very long on GroupOnFilter!
-        {
-          SMESH::smIdType_array_var nbByType = mesh_i->GetNbElementsByType();
-          isEmpty = ( nbByType[ grp->GetType() ] == 0 );
+        SMESH::SMESH_GroupBase_var grp = SMESH::SMESH_GroupBase::_narrow( obj );
+        if ( !grp->_is_nil() ) {
+          UpdateGroupIcon(grp);
+          continue;
         }
         else
         {
@@ -1013,17 +1010,46 @@ void SMESH_Gen_i::UpdateIcons( SMESH::SMESH_Mesh_ptr theMesh )
 
       if ( isEmpty )
         SetPixMap( so, "ICON_SMESH_TREE_MESH_WARN");
-      else if ( !isGroup )
-        SetPixMap( so, "ICON_SMESH_TREE_MESH" );
-      else if ( isGroupOnFilter )
-        SetPixMap( so, "ICON_SMESH_TREE_GROUP_ON_FILTER" );
       else
-        SetPixMap( so, "ICON_SMESH_TREE_GROUP" );
+        SetPixMap( so, "ICON_SMESH_TREE_MESH" );
 
     } // loop on sub-meshes or groups
   } // loop on roots
 }
 
+//=======================================================================
+//function : UpdateGroupIcon
+//purpose  : update icon of a group
+//=======================================================================
+
+void SMESH_Gen_i::UpdateGroupIcon( SMESH::SMESH_GroupBase_ptr theGroup )
+{
+  SALOMEDS::SObject_wrap so = ObjectToSObject( theGroup );
+  if ( so->_is_nil() )
+    return;
+
+  SMESH::SMESH_GroupOnFilter_var gof = SMESH::SMESH_GroupOnFilter::_narrow( theGroup );
+  const bool isGroupOnFilter = !gof->_is_nil();
+
+  bool isEmpty = false;
+  if ( isGroupOnFilter ) // GetTypes() can be very long on GroupOnFilter!
+  {
+    SMESH::smIdType_array_var nbByType = theGroup->GetMesh()->GetNbElementsByType();
+    isEmpty = ( nbByType[ theGroup->GetType() ] == 0 );
+  }
+  else
+  {
+    isEmpty = ( theGroup->Size() == 0 );
+  }
+
+  if ( isEmpty )
+    SetPixMap( so, "ICON_SMESH_TREE_MESH_WARN");
+  else if ( isGroupOnFilter )
+    SetPixMap( so, "ICON_SMESH_TREE_GROUP_ON_FILTER" );
+  else
+    SetPixMap( so, "ICON_SMESH_TREE_GROUP" );
+}
+
 //=======================================================================
 //function : HighLightInvalid
 //purpose  : change font color of a object in the Object Browser

src/SMESH_I/SMESH_Group_i.cxx

@@ -261,16 +261,19 @@ void SMESH_Group_i::Clear()
   if ( myPreMeshInfo )
     myPreMeshInfo->FullLoadFromFile();
 
+  SMESH::SMESH_Group_var me = _this();
+
   // Update Python script
-  SMESH::TPythonDump() << SMESH::SMESH_Group_var(_this()) << ".Clear()";
+  SMESH::TPythonDump() << me << ".Clear()";
 
   // Clear the group
   SMESHDS_Group* aGroupDS = dynamic_cast<SMESHDS_Group*>( GetGroupDS() );
   if (aGroupDS) {
     aGroupDS->Clear();
-    return;
+    SMESH_Gen_i* aGen = GetMeshServant()->GetGen();
+    aGen->UpdateGroupIcon(me);
+    Modified(); // notify dependent Filter with FT_BelongToMeshGroup criterion
   }
-  Modified(); // notify dependent Filter with FT_BelongToMeshGroup criterion
 }
 
 //=============================================================================
@@ -301,20 +304,28 @@ SMESH::smIdType  SMESH_Group_i::Add( const SMESH::smIdType_array& theIDs )
   if ( myPreMeshInfo )
     myPreMeshInfo->FullLoadFromFile();
 
+  SMESH::SMESH_Group_var me = _this();
+
   // Update Python script
-  SMESH::TPythonDump() << "nbAdd = " << SMESH::SMESH_Group_var(_this()) << ".Add( " << theIDs << " )";
+  SMESH::TPythonDump() << "nbAdd = " << me << ".Add( " << theIDs << " )";
 
   // Add elements to the group
   SMESHDS_Group* aGroupDS = dynamic_cast<SMESHDS_Group*>( GetGroupDS() );
   if (aGroupDS) {
+    bool wasEmpty = aGroupDS->IsEmpty();
     int nbAdd = 0;
     for ( CORBA::ULong i = 0; i < theIDs.length(); i++) {
       int anID = (int) theIDs[i];
       if ( aGroupDS->Add( anID ))
         nbAdd++;
     }
-    if ( nbAdd )
+    if ( nbAdd ) {
       Modified(); // notify dependent Filter with FT_BelongToMeshGroup criterion
+      if (wasEmpty) {
+        SMESH_Gen_i* aGen = GetMeshServant()->GetGen();
+        aGen->UpdateGroupIcon(me);
+      }
+    }
     return nbAdd;
   }
   MESSAGE("attempt to add elements to a vague group");
@@ -332,9 +343,10 @@ SMESH::smIdType  SMESH_Group_i::Remove( const SMESH::smIdType_array& theIDs )
   if ( myPreMeshInfo )
     myPreMeshInfo->FullLoadFromFile();
 
+  SMESH::SMESH_Group_var me = _this();
+
   // Update Python script
-  SMESH::TPythonDump() << "nbDel = " << SMESH::SMESH_Group_var(_this())
-                << ".Remove( " << theIDs << " )";
+  SMESH::TPythonDump() << "nbDel = " << me << ".Remove( " << theIDs << " )";
 
   // Remove elements from the group
   SMESHDS_Group* aGroupDS = dynamic_cast<SMESHDS_Group*>( GetGroupDS() );
@@ -345,8 +357,13 @@ SMESH::smIdType  SMESH_Group_i::Remove( const SMESH::smIdType_array& theIDs )
       if ( aGroupDS->Remove( anID ))
         nbDel++;
     }
-    if ( nbDel )
+    if ( nbDel ) {
       Modified(); // notify dependent Filter with FT_BelongToMeshGroup criterion
+      if (aGroupDS->IsEmpty()) {
+        SMESH_Gen_i* aGen = GetMeshServant()->GetGen();
+        aGen->UpdateGroupIcon(me);
+      }
+    }
     return nbDel;
   }
   MESSAGE("attempt to remove elements from a vague group");
@@ -361,24 +378,27 @@ SMESH::smIdType  SMESH_Group_i::Remove( const SMESH::smIdType_array& theIDs )
 
 typedef bool (SMESHDS_Group::*TFunChangeGroup)(const smIdType);
 
-CORBA::Long 
-ChangeByPredicate( SMESH::Predicate_i*       thePredicate,
-                   SMESHDS_GroupBase*        theGroupBase,
-                   SMESH::NotifyerAndWaiter* theGroupImpl,
-                   TFunChangeGroup           theFun)
+CORBA::Long ChangeByPredicate( SMESH::Predicate_i* thePredicate,
+                               SMESHDS_GroupBase*  theGroupBase,
+                               SMESH_GroupBase_i*  theGroup,
+                               TFunChangeGroup     theFun)
 {
   CORBA::Long aNb = 0;
-  if(SMESHDS_Group* aGroupDS = dynamic_cast<SMESHDS_Group*>(theGroupBase)){
+  if (SMESHDS_Group* aGroupDS = dynamic_cast<SMESHDS_Group*>(theGroupBase)) {
     SMESH::Controls::Filter::TIdSequence aSequence;
     const SMDS_Mesh* aMesh = theGroupBase->GetMesh();
     SMESH::Filter_i::GetElementsId(thePredicate,aMesh,aSequence);
     
     CORBA::Long i = 0, iEnd = aSequence.size();
-    for(; i < iEnd; i++)
-      if((aGroupDS->*theFun)(aSequence[i]))
+    for (; i < iEnd; i++)
+      if ((aGroupDS->*theFun)(aSequence[i]))
         aNb++;
-    if ( aNb )
-      theGroupImpl->Modified();
+    if ( aNb ) {
+      theGroup->Modified();
+      SMESH_Gen_i* aGen = theGroup->GetMeshServant()->GetGen();
+      SMESH::SMESH_GroupBase_var aGroup = theGroup->_this();
+      aGen->UpdateGroupIcon(aGroup);
+    }
     return aNb;
   }
   return aNb;
@@ -428,11 +448,16 @@ SMESH::smIdType  SMESH_Group_i::AddFrom( SMESH::SMESH_IDSource_ptr theSource )
         aGroupDS->SMDSGroup().Add( elemIt->next() );
   }
 
-  // Update Python script
-  pd << "nbAdd = " << SMESH::SMESH_Group_var(_this()) << ".AddFrom( " << theSource << " )";
+  SMESH::SMESH_Group_var me = _this();
 
-  if ( prevNb != Size() )
+  // Update Python script
+  pd << "nbAdd = " << me << ".AddFrom( " << theSource << " )";
+
+  if ( prevNb != Size() ) {
+    SMESH_Gen_i* aGen = GetMeshServant()->GetGen();
+    aGen->UpdateGroupIcon(me);
     Modified(); // notify dependent Filter with FT_BelongToMeshGroup criterion
+  }
 
   return Size() - prevNb;
 }
@@ -830,17 +855,18 @@ void SMESH_GroupOnFilter_i::SetFilter(SMESH::Filter_ptr theFilter)
     myFilter->SetMesh( SMESH::SMESH_Mesh::_nil() ); // to UnRegister() the mesh
   }
 
+  SMESH_Gen_i* aGen = GetMeshServant()->GetGen();
+  SMESH::SMESH_GroupOnFilter_var me = _this();
+
   if ( SMESHDS_GroupOnFilter* grDS = dynamic_cast< SMESHDS_GroupOnFilter*>( GetGroupDS() ))
   {
     grDS->SetPredicate( GetPredicate( myFilter ));
+    aGen->UpdateGroupIcon(me);
     Modified(); // notify dependent Filter with FT_BelongToMeshGroup criterion
   }
 
-  SMESH::SMESH_GroupOnFilter_var me = _this();
-
   // mark the group valid after edition
-  GetMeshServant()->GetGen()->HighLightInvalid( me, false );
-
+  aGen->HighLightInvalid( me, false );
 
   SMESH::TPythonDump()<< me <<".SetFilter( "<< theFilter <<" )";
 }
@@ -1078,4 +1104,3 @@ void SMESH_GroupOnFilter_i::OnBaseObjModified(NotifyerAndWaiter* /*filter*/, boo
   if ( SMESHDS_GroupOnFilter* grDS = dynamic_cast< SMESHDS_GroupOnFilter*>( GetGroupDS() ))
     grDS->SetPredicate( GetPredicate( myFilter )); // group resets its cache
 }
-

src/SMESH_I/SMESH_Mesh_i.cxx

@@ -1414,6 +1414,9 @@ SMESH::SMESH_Group_ptr SMESH_Mesh_i::IntersectGroups( SMESH::SMESH_GroupBase_ptr
         resGroupDS->SMDSGroup().Add( e );
     }
   }
+
+  GetGen()->UpdateGroupIcon(aResGrp);
+
   // Update Python script
   pyDump << aResGrp << " = " << SMESH::SMESH_Mesh_var(_this()) << ".IntersectGroups( "
          << theGroup1 << ", " << theGroup2 << ", '" << theName << "')";
@@ -1500,6 +1503,8 @@ SMESH_Mesh_i::IntersectListOfGroups(const SMESH::ListOfGroups& theGroups,
       resGroupDS->SMDSGroup().Add( e );
   }
 
+  GetGen()->UpdateGroupIcon(aResGrp);
+
   // Update Python script
   pyDump << aResGrp << " = " << SMESH::SMESH_Mesh_var( _this() )
          << ".IntersectListOfGroups( " << theGroups << ", '" << theName << "' )";
@@ -1561,6 +1566,9 @@ SMESH::SMESH_Group_ptr SMESH_Mesh_i::CutGroups( SMESH::SMESH_GroupBase_ptr theGr
         resGroupDS->SMDSGroup().Add( e );
     }
   }
+
+  GetGen()->UpdateGroupIcon(aResGrp);
+
   // Update Python script
   pyDump << aResGrp << " = " << SMESH::SMESH_Mesh_var(_this()) << ".CutGroups( "
          << theGroup1 << ", " << theGroup2 << ", '" << theName << "')";
@@ -1661,6 +1669,8 @@ SMESH_Mesh_i::CutListOfGroups(const SMESH::ListOfGroups& theMainGroups,
       resGroupDS->SMDSGroup().Add( e );
   }
 
+  GetGen()->UpdateGroupIcon(aResGrp);
+
   // Update Python script
   pyDump << aResGrp << " = " << SMESH::SMESH_Mesh_var( _this() )
          << ".CutListOfGroups( " << theMainGroups << ", "
@@ -3473,7 +3483,6 @@ void SMESH_Mesh_i::onHypothesisModified(int theHypID, bool theUpdateIcons)
 
 void SMESH_Mesh_i::SetImpl(::SMESH_Mesh * impl)
 {
-  MESSAGE("SMESH_Mesh_i::SetImpl");
   _impl = impl;
   if ( _impl )
     _impl->SetCallUp( new TCallUp_i(this));
@@ -3487,7 +3496,6 @@ void SMESH_Mesh_i::SetImpl(::SMESH_Mesh * impl)
 
 ::SMESH_Mesh & SMESH_Mesh_i::GetImpl()
 {
-  MESSAGE("SMESH_Mesh_i::GetImpl()");
   return *_impl;
 }
 
@@ -4576,7 +4584,7 @@ void SMESH_Mesh_i::ExportStructuredCGNS( SMESH::SMESH_IDSource_ptr meshPart,
     meshName = name.in();
   }
   SMESH_TRY;
-  
+
   SMESH::SMESH_Mesh_var mesh = meshPart->GetMesh();
   SMESH_Mesh_i*       mesh_i = SMESH::DownCast<SMESH_Mesh_i*>( mesh );
   mesh_i->Load();

src/SMESH_I/SMESH_Meshio.cxx

@@ -365,7 +365,9 @@ QString SMESH_Meshio::GetConvertOptArgs() const
     return "-o ansys";
   else if (mySelectedFilter.startsWith("Gmsh 2"))
     return "-o gmsh22";
-  else if (mySelectedFilter.startsWith("Gmsh 4"))
+  else if (mySelectedFilter.startsWith("Gmsh 4.0"))
+    return "-o gmsh40";
+  else if (mySelectedFilter.startsWith("Gmsh 4.1"))
     return "-o gmsh";
 
   return {};

src/SMESH_SWIG/smeshBuilder.py

@@ -675,6 +675,14 @@ class smeshBuilder( SMESH._objref_SMESH_Gen, object ):
         global notebook
         notebook = salome_notebook.NoteBook( theIsEnablePublish )
 
+    def ReloadMeshFromFile(self, theMesh):
+        """
+        Desc for method,
+        """
+
+        aSmeshMesh = SMESH._objref_SMESH_Gen.ReloadMeshFromFile(self, theMesh)
+        aMesh = Mesh(self, self.geompyD, aSmeshMesh)
+        return aMesh
 
     def CreateMeshesFromUNV( self,theFileName ):
         """

src/StdMeshers.test/CMakeLists.txt

@@ -0,0 +1,70 @@
+# Copyright (C) 2012-2024  CEA, EDF, OPEN CASCADE
+#
+# This library is free software; you can redistribute it and/or
+# modify it under the terms of the GNU Lesser General Public
+# License as published by the Free Software Foundation; either
+# version 2.1 of the License, or (at your option) any later version.
+#
+# This library is distributed in the hope that it will be useful,
+# but WITHOUT ANY WARRANTY; without even the implied warranty of
+# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+# Lesser General Public License for more details.
+#
+# You should have received a copy of the GNU Lesser General Public
+# License along with this library; if not, write to the Free Software
+# Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
+#
+# See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
+#
+
+INCLUDE(tests.set)
+
+SALOME_GENERATE_TESTS_ENVIRONMENT(_test_env)
+
+SET(TEST_INSTALL_DIRECTORY ${SMESH_TEST_DIR}/other)
+
+INCLUDE_DIRECTORIES(
+  ${OpenCASCADE_INCLUDE_DIR}
+  ${Boost_INCLUDE_DIRS}
+  ${MEDCOUPLING_INCLUDE_DIRS}
+  ${CPPUNIT_INCLUDE_DIRS}
+  ${PROJECT_SOURCE_DIR}/src/StdMeshers
+  ${PROJECT_SOURCE_DIR}/src/SMESHUtils
+  ${PROJECT_SOURCE_DIR}/src/SMESH
+  ${PROJECT_SOURCE_DIR}/src/SMESHDS
+  ${PROJECT_SOURCE_DIR}/src/SMDS
+  ${PROJECT_SOURCE_DIR}/src/Controls
+)
+# additional preprocessor / compiler flags
+ADD_DEFINITIONS(
+  ${OpenCASCADE_DEFINITIONS}
+  ${BOOST_DEFINITIONS}
+  ${CPPUNIT_DEFINITIONS}
+)
+
+IF(SALOME_SMESH_USE_TBB)
+  SET(TBB_LIBS ${TBB_LIBRARIES})
+ENDIF(SALOME_SMESH_USE_TBB)
+
+FOREACH(_test ${UNIT_TESTS})
+  GET_FILENAME_COMPONENT(testname ${_test} NAME_WE)
+  SET(testname "TESTS_${testname}")
+
+  add_executable(${_test} ${_test}.cxx)
+  target_link_libraries(${_test} StdMeshers ${CPPUNIT_LIBRARIES} )
+
+  ADD_TEST(NAME ${testname}
+           COMMAND ${CMAKE_CURRENT_SOURCE_DIR}/${_test} )
+  SET_TESTS_PROPERTIES(${testname} PROPERTIES ENVIRONMENT "${tests_env}" LABELS "tests")
+ENDFOREACH()
+
+IF(WIN32)
+  FOREACH(_test ${UNIT_TESTS})
+    INSTALL(FILES ${CMAKE_CURRENT_BINARY_DIR}/${CMAKE_BUILD_TYPE}/${_test}${CMAKE_EXECUTABLE_SUFFIX} PERMISSIONS OWNER_EXECUTE OWNER_WRITE OWNER_READ DESTINATION ${TEST_INSTALL_DIRECTORY})
+  ENDFOREACH()
+ELSE()
+  FOREACH(_test ${UNIT_TESTS})
+    INSTALL(FILES ${CMAKE_CURRENT_BINARY_DIR}/${_test} PERMISSIONS OWNER_EXECUTE OWNER_WRITE OWNER_READ DESTINATION ${TEST_INSTALL_DIRECTORY})
+  ENDFOREACH()
+ENDIF(WIN32)
+

src/StdMeshers.test/HexahedronCanonicalShapesTest.cxx

@@ -0,0 +1,245 @@
+// Copyright (C) 2016-2024  CEA, EDF
+//
+// This library is free software; you can redistribute it and/or
+// modify it under the terms of the GNU Lesser General Public
+// License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+//
+// This library is distributed in the hope that it will be useful,
+// but WITHOUT ANY WARRANTY; without even the implied warranty of
+// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+// Lesser General Public License for more details.
+//
+// You should have received a copy of the GNU Lesser General Public
+// License along with this library; if not, write to the Free Software
+// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
+//
+// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
+//
+//  File   : HexahedronCanonicalShapesTest.cxx
+//  Module : SMESH
+//  Purpose: Implement unit tests for StdMeshers_Cartesian_3D_Hexahedron class to reproduce bugs that manifest in integration tests.
+//            The main difference between this unit test and integration tests is the fine grained control we have over the class methods and the hability to diagnose/solve bugs before the code goes into production enviroment.
+//            This test class can be used as reference for the development of future tests in other stdMesh algorithms
+
+#include "StdMeshers_Cartesian_3D_Hexahedron.hxx"
+#include "StdMeshers_CartesianParameters3D.hxx"
+
+// CPP TEST
+#include <cppunit/TestAssert.h>
+
+// OCC
+#include <BRep_Builder.hxx>
+#include <BRepTools.hxx>
+#include <BRepPrimAPI_MakeBox.hxx>
+#include <BRepPrimAPI_MakeCylinder.hxx>
+#include <BRepPrimAPI_MakeSphere.hxx>
+#include <BRepPrimAPI_MakeCone.hxx>
+
+#include <iostream>
+#include <memory>
+
+using namespace StdMeshers::Cartesian3D;
+
+// Helper functions!
+// Build Grid
+//      Require building mesh
+//      Require building shape.
+
+/*!
+  * \brief Mock mesh
+  */
+struct SMESH_Mesh_Test: public SMESH_Mesh
+{
+  SMESH_Mesh_Test() {
+    _isShapeToMesh = (_id = 0);
+    _meshDS  = new SMESHDS_Mesh( _id, true );
+  }
+};
+
+/*!
+  * \brief Mock Hypothesis
+  */
+struct CartesianHypo: public StdMeshers_CartesianParameters3D
+{
+  CartesianHypo() : StdMeshers_CartesianParameters3D(0/*zero hypoId*/, nullptr/*NULL generator*/)
+  {
+  }
+};
+
+/*!
+  * \brief Initialize the grid and intesersectors of grid with the geometry
+  */
+void GridInitAndIntersectWithShape (Grid& grid,
+                                    double gridSpacing,
+                                    double theSizeThreshold,
+                                    const TopoDS_Shape theShape,
+                                    TEdge2faceIDsMap& edge2faceIDsMap,
+                                    const int theNumOfThreads)
+{
+  // Canonical axes(i,j,k)
+  double axisDirs[9] = {1.,0.,0.,0.,1.,0.,0.,0.,1.};
+  std::vector<std::string> grdSpace = { std::to_string(gridSpacing) };
+  std::vector<double> intPnts;
+
+  std::unique_ptr<CartesianHypo> aHypo ( new CartesianHypo() );
+  aHypo->SetGridSpacing(grdSpace, intPnts, 0 ); // Spacing in dir 0
+  aHypo->SetGridSpacing(grdSpace, intPnts, 1 ); // Spacing in dir 1
+  aHypo->SetGridSpacing(grdSpace, intPnts, 2 ); // Spacing in dir 2
+  aHypo->SetSizeThreshold(theSizeThreshold);    // set threshold
+  aHypo->SetAxisDirs(axisDirs);
+
+  grid.GridInitAndInterserctWithShape(theShape, edge2faceIDsMap, aHypo.get(), theNumOfThreads, false);
+}
+
+/*!
+  * \brief Test runner
+  */
+bool testShape (const TopoDS_Shape theShape,
+                const bool toAddEdges,
+                const bool toCreateFaces,
+                const double theGridSpacing,
+                const double theSizeThreshold,
+                const int theNbCreatedExpected)
+{
+  std::unique_ptr<SMESH_Mesh> aMesh( new SMESH_Mesh_Test() );
+  aMesh->ShapeToMesh( theShape );
+  SMESH_MesherHelper helper( *aMesh );
+
+  Grid grid;
+  grid._helper = &helper;
+  grid._toAddEdges = toAddEdges;
+  grid._toCreateFaces = toCreateFaces;
+  grid._toConsiderInternalFaces = false;
+  grid._toUseThresholdForInternalFaces = false;
+  grid._toUseQuanta = false;
+  grid._sizeThreshold = theSizeThreshold;
+
+  TEdge2faceIDsMap edge2faceIDsMap;
+  GridInitAndIntersectWithShape( grid, theGridSpacing, theSizeThreshold,
+                                 theShape, edge2faceIDsMap, 1 );
+
+  SMESH_subMesh * aSubMesh = aMesh->GetSubMesh(theShape);
+  aSubMesh->DependsOn(); // init sub-meshes
+
+  Hexahedron hex( &grid );
+  int nbAdded = hex.MakeElements( helper, edge2faceIDsMap, 1 );
+  if (nbAdded != theNbCreatedExpected) {
+    std::stringstream buffer;
+    buffer << "Number of computed elements does not match: obtained " << nbAdded << " != expected " << theNbCreatedExpected;
+    //CPPUNIT_ASSERT_MESSAGE(buffer.str().c_str(), nbAdded == theNbCreatedExpected );
+    //MESSAGE(buffer.str().c_str());
+    //CppUnitTestFramework::Logger::WriteMessage(buffer.str().c_str());
+    std::cerr << buffer.str() << std::endl;
+    return false;
+  }
+
+  return true;
+}
+
+/*!
+  * \brief Test some primitive shapes
+  */
+bool testPrimitives()
+{
+  bool isOK = true;
+
+  // Test fitting of a box
+  BRepPrimAPI_MakeBox aMakeBox (10, 20, 30);
+  aMakeBox.Build();
+  CPPUNIT_ASSERT_MESSAGE( "Could not create the box!", aMakeBox.IsDone() );
+  TopoDS_Shape aShape = aMakeBox.Shape();
+
+  // Test exact fitting of a box
+  if (!testShape (aShape, /*toAddEdges*/false, /*toCreateFaces*/false,
+                  /*gridSpacing*/10, /*theSizeThreshold*/4, /*theNbCreatedExpected*/6))
+    isOK = false;
+  if (!testShape (aShape, /*toAddEdges*/true, /*toCreateFaces*/false,
+                  /*gridSpacing*/10, /*theSizeThreshold*/4, /*theNbCreatedExpected*/6))
+    isOK = false;
+
+  if (!testShape (aShape, /*toAddEdges*/false, /*toCreateFaces*/true,
+                  /*gridSpacing*/10, /*theSizeThreshold*/4, /*theNbCreatedExpected*/6))
+    isOK = false;
+
+  if (!testShape (aShape, /*toAddEdges*/false, /*toCreateFaces*/false,
+                  /*gridSpacing*/5, /*theSizeThreshold*/4, /*theNbCreatedExpected*/48))
+    isOK = false;
+
+  // Test not exact fitting of a box
+  if (!testShape (aShape, /*toAddEdges*/false, /*toCreateFaces*/false,
+                  /*gridSpacing*/7, /*theSizeThreshold*/4, /*theNbCreatedExpected*/12))
+    isOK = false;
+
+  // Test fitting of a cylinder
+  gp_Ax2 anAxes (gp::Origin(), gp::DZ());
+  BRepPrimAPI_MakeCylinder aMakeCyl (anAxes, 20., 30.);
+  aMakeCyl.Build();
+  CPPUNIT_ASSERT_MESSAGE( "Could not create the cylinder!", aMakeCyl.IsDone() );
+  aShape = aMakeCyl.Shape();
+
+  // test for different threshold values
+  if (!testShape (aShape, /*toAddEdges*/false, /*toCreateFaces*/false,
+                  /*gridSpacing*/10, /*theSizeThreshold*/4, /*theNbCreatedExpected*/48))
+    isOK = false;
+  if (!testShape (aShape, /*toAddEdges*/false, /*toCreateFaces*/false,
+                  /*gridSpacing*/10, /*theSizeThreshold*/2, /*theNbCreatedExpected*/36))
+    isOK = false;
+
+  // Test fitting of a sphere
+  BRepPrimAPI_MakeSphere aMakeSph (anAxes, 30.);
+  aMakeSph.Build();
+  CPPUNIT_ASSERT_MESSAGE( "Could not create the sphere!", aMakeSph.IsDone() );
+  aShape = aMakeSph.Shape();
+
+  // test for different threshold values
+  if (!testShape (aShape, /*toAddEdges*/false, /*toCreateFaces*/false,
+                  /*gridSpacing*/10, /*theSizeThreshold*/4, /*theNbCreatedExpected*/136))
+    isOK = false;
+  if (!testShape (aShape, /*toAddEdges*/false, /*toCreateFaces*/false,
+                  /*gridSpacing*/10, /*theSizeThreshold*/2, /*theNbCreatedExpected*/88))
+    isOK = false;
+
+  // Test fitting of a cone
+  BRepPrimAPI_MakeCone aMakeCon (anAxes, 30., 0., 40.);
+  aMakeCon.Build();
+  CPPUNIT_ASSERT_MESSAGE( "Could not create the cone!", aMakeCon.IsDone() );
+  aShape = aMakeCon.Shape();
+
+  // test for different threshold values
+  if (!testShape (aShape, /*toAddEdges*/false, /*toCreateFaces*/false,
+                  /*gridSpacing*/10, /*theSizeThreshold*/100, /*theNbCreatedExpected*/72))
+    isOK = false;
+  if (!testShape (aShape, /*toAddEdges*/false, /*toCreateFaces*/false,
+                  /*gridSpacing*/10, /*theSizeThreshold*/4, /*theNbCreatedExpected*/40))
+    isOK = false;
+  if (!testShape (aShape, /*toAddEdges*/false, /*toCreateFaces*/false,
+                  /*gridSpacing*/10, /*theSizeThreshold*/1.5, /*theNbCreatedExpected*/32))
+    isOK = false;
+
+  // truncated cone
+  aMakeCon = BRepPrimAPI_MakeCone(anAxes, 30., 15., 20.);
+  aMakeCon.Build();
+  CPPUNIT_ASSERT_MESSAGE( "Could not create the cone!", aMakeCon.IsDone() );
+  aShape = aMakeCon.Shape();
+
+  // test for different threshold values
+  if (!testShape (aShape, /*toAddEdges*/false, /*toCreateFaces*/false,
+                  /*gridSpacing*/10, /*theSizeThreshold*/100, /*theNbCreatedExpected*/56))
+    isOK = false;
+  if (!testShape (aShape, /*toAddEdges*/false, /*toCreateFaces*/false,
+                  /*gridSpacing*/10, /*theSizeThreshold*/4, /*theNbCreatedExpected*/36))
+    isOK = false;
+  if (!testShape (aShape, /*toAddEdges*/false, /*toCreateFaces*/false,
+                  /*gridSpacing*/10, /*theSizeThreshold*/1.5, /*theNbCreatedExpected*/28))
+    isOK = false;
+
+  return isOK;
+}
+
+// Entry point for test
+int main()
+{
+  bool isOK = testPrimitives();
+  return isOK ? 0 : 1;
+}

src/StdMeshers.test/HexahedronTest.cxx

@@ -0,0 +1,204 @@
+// Copyright (C) 2016-2024  CEA, EDF
+//
+// This library is free software; you can redistribute it and/or
+// modify it under the terms of the GNU Lesser General Public
+// License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+//
+// This library is distributed in the hope that it will be useful,
+// but WITHOUT ANY WARRANTY; without even the implied warranty of
+// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+// Lesser General Public License for more details.
+//
+// You should have received a copy of the GNU Lesser General Public
+// License along with this library; if not, write to the Free Software
+// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
+//
+// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
+//
+//  File   : HexahedronTest.cxx
+//  Module : SMESH
+//  Purpose: Implement unit tests for StdMeshers_Cartesian_3D_Hexahedron class to reproduce bugs that manifest in integration tests.
+//            The main difference between this unit test and integration tests is the fine grained control we have over the class methods and the hability to diagnose/solve bugs before the code goes into production enviroment.
+//            This test class can be used as reference for the development of future tests in other stdMesh algorithms
+
+#include "StdMeshers_Cartesian_3D_Hexahedron.hxx"
+#include "StdMeshers_CartesianParameters3D.hxx"
+
+// CPP TEST
+#include <cppunit/TestAssert.h>
+
+// OCC
+#include <BRep_Builder.hxx>
+#include <BRepTools.hxx>
+#include <BRepPrimAPI_MakeBox.hxx>
+#include <BRepPrimAPI_MakeCylinder.hxx>
+#include <BRepPrimAPI_MakeSphere.hxx>
+#include <BRepPrimAPI_MakeCone.hxx>
+
+#include <iostream>
+#include <memory>
+#include <numeric>
+
+using namespace StdMeshers::Cartesian3D;
+
+// Helper functions!
+// Build Grid
+//      Require building mesh
+//      Require building shape.
+
+/*!
+  * \brief Mock mesh
+  */
+struct SMESH_Mesh_Test: public SMESH_Mesh
+{
+  SMESH_Mesh_Test() {
+    _isShapeToMesh = (_id = 0);
+    _meshDS  = new SMESHDS_Mesh( _id, true );
+  }
+};
+
+/*!
+  * \brief Mock Hypothesis
+  */
+struct CartesianHypo: public StdMeshers_CartesianParameters3D
+{
+  CartesianHypo() : StdMeshers_CartesianParameters3D(0/*zero hypoId*/, nullptr/*NULL generator*/)
+  {
+  }
+};
+
+/*!
+  * \brief Shape loader
+  */
+void loadBrepShape( std::string shapeName, TopoDS_Shape & shape )
+{
+  BRep_Builder b;
+  BRepTools::Read(shape, shapeName.c_str(), b);
+}
+
+/*!
+  * \brief Initialize the grid and intesersectors of grid with the geometry
+  */
+void GridInitAndIntersectWithShape (Grid& grid,
+                                    double gridSpacing,
+                                    double theSizeThreshold,
+                                    const TopoDS_Shape theShape,
+                                    TEdge2faceIDsMap& edge2faceIDsMap,
+                                    const int theNumOfThreads)
+{
+  // Canonical axes(i,j,k)
+  double axisDirs[9] = {1.,0.,0., 0.,1.,0., 0.,0.,1.};
+  std::vector<std::string> grdSpace = { std::to_string(gridSpacing) };
+  std::vector<double> intPnts;
+
+  std::unique_ptr<CartesianHypo> aHypo ( new CartesianHypo() );
+  aHypo->SetAxisDirs(axisDirs);
+  aHypo->SetGridSpacing(grdSpace, intPnts, 0 ); // Spacing in dir 0
+  aHypo->SetGridSpacing(grdSpace, intPnts, 1 ); // Spacing in dir 1
+  aHypo->SetGridSpacing(grdSpace, intPnts, 2 ); // Spacing in dir 2
+  aHypo->SetSizeThreshold(theSizeThreshold);    // set threshold
+
+  grid.GridInitAndInterserctWithShape(theShape, edge2faceIDsMap, aHypo.get(), theNumOfThreads, false);
+}
+
+/*!
+  * \brief Reproduce conditions of TBPERF_GRIDS_PERF_SMESH_M1 test to detect and solve segfault in unit test.
+  */
+bool testNRTM1()
+{
+  TopoDS_Shape aShape;
+  loadBrepShape( "data/HexahedronTest/NRTM1.brep", aShape );
+  CPPUNIT_ASSERT_MESSAGE( "Could not load the brep shape!", !aShape.IsNull() );
+
+  const auto numOfCores = std::thread::hardware_concurrency() == 0 ? 1 : std::thread::hardware_concurrency();
+  std::vector<int> numberOfThreads(numOfCores);
+  std::iota (std::begin(numberOfThreads), std::end(numberOfThreads), 1);
+
+  for (auto nThreads : numberOfThreads )
+  {
+    for (size_t i = 0; i < 10 /*trials*/; i++)
+    {
+      std::unique_ptr<SMESH_Mesh> aMesh( new SMESH_Mesh_Test() );
+      aMesh->ShapeToMesh( aShape );
+      SMESH_MesherHelper helper( *aMesh );
+
+      Grid grid;
+      grid._helper = &helper;
+      grid._toAddEdges = false;
+      grid._toCreateFaces = false;
+      grid._toConsiderInternalFaces = false;
+      grid._toUseThresholdForInternalFaces = false;
+      grid._toUseQuanta = false;
+      grid._sizeThreshold = 4.0;
+
+      TEdge2faceIDsMap edge2faceIDsMap;
+      GridInitAndIntersectWithShape( grid, 1.0, 4.0, aShape, edge2faceIDsMap, nThreads );
+
+      SMESH_subMesh * aSubMesh = aMesh->GetSubMesh(aShape);
+      aSubMesh->DependsOn(); // init sub-meshes
+
+      Hexahedron hex( &grid );
+      int nbAdded = hex.MakeElements( helper, edge2faceIDsMap, nThreads );
+      CPPUNIT_ASSERT_MESSAGE( "Number of computed elements does not match", nbAdded == 1024 );
+    }
+  }
+  return true;
+}
+
+/*!
+  * \brief Reproduce conditions of TBPERF_GRIDS_PERF_SMESH_J4 test to detect and solve segfault in unit test.
+  */
+bool testNRTJ4()
+{
+  TopoDS_Shape aShape;
+  loadBrepShape( "data/HexahedronTest/NRTMJ4.brep", aShape );
+  CPPUNIT_ASSERT_MESSAGE( "Could not load the brep shape!", !aShape.IsNull() );
+
+  const auto numOfCores = std::thread::hardware_concurrency() == 0 ? 1 : std::thread::hardware_concurrency()/2;
+  std::vector<int> numberOfThreads(numOfCores);
+  std::iota (std::begin(numberOfThreads), std::end(numberOfThreads), 1);
+
+  // Test with face creation
+  for (auto nThreads : numberOfThreads )
+  {
+    for (size_t i = 0; i < 10 /*trials*/; i++)
+    {
+      std::unique_ptr<SMESH_Mesh> aMesh( new SMESH_Mesh_Test() );
+      aMesh->ShapeToMesh( aShape );
+      SMESH_MesherHelper helper( *aMesh );
+
+      Grid grid;
+      grid._helper = &helper;
+      grid._toAddEdges = false;
+      grid._toConsiderInternalFaces = false;
+      grid._toUseThresholdForInternalFaces = false;
+      grid._toUseQuanta = false;
+      double testThreshold = 1.000001;
+      grid._toCreateFaces = true;
+      grid._sizeThreshold = testThreshold;
+
+      TEdge2faceIDsMap edge2faceIDsMap;
+      GridInitAndIntersectWithShape( grid, 2.0, testThreshold,
+                                     aShape, edge2faceIDsMap, nThreads );
+
+      SMESH_subMesh * aSubMesh = aMesh->GetSubMesh(aShape);
+      aSubMesh->DependsOn(); // init sub-meshes
+
+      Hexahedron hex( &grid );
+      int nbAdded = hex.MakeElements( helper, edge2faceIDsMap, nThreads );
+      CPPUNIT_ASSERT_MESSAGE( "Number of computed elements does not match", nbAdded == 35150 );
+    }
+  }
+  return true;
+}
+
+// Entry point for test
+int main()
+{
+  bool isOK = testNRTM1();
+  if (!testNRTJ4())
+    isOK = false;
+
+  return isOK ? 0 : 1;
+}

src/StdMeshers.test/tests.set

@@ -0,0 +1,28 @@
+# Copyright (C) 2015-2024  CEA, EDF, OPEN CASCADE
+#
+# This library is free software; you can redistribute it and/or
+# modify it under the terms of the GNU Lesser General Public
+# License as published by the Free Software Foundation; either
+# version 2.1 of the License, or (at your option) any later version.
+#
+# This library is distributed in the hope that it will be useful,
+# but WITHOUT ANY WARRANTY; without even the implied warranty of
+# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+# Lesser General Public License for more details.
+#
+# You should have received a copy of the GNU Lesser General Public
+# License along with this library; if not, write to the Free Software
+# Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
+#
+# See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
+#
+
+# The following tests cannot be executed with 'make test' because they use
+# external meshing plug-ins.
+# On the other hand these tests can be executed with 'salome test'.
+# ---------------------------------------------------------------------------
+
+SET(UNIT_TESTS
+  HexahedronTest
+  HexahedronCanonicalShapesTest
+  )

src/StdMeshers/CMakeLists.txt

@@ -119,6 +119,8 @@ SET(StdMeshers_HEADERS
   StdMeshers_ViscousLayers2D.hxx
   StdMeshers_Projection_1D2D.hxx
   StdMeshers_CartesianParameters3D.hxx
+  StdMeshers_Cartesian_3D_Grid.hxx
+  StdMeshers_Cartesian_3D_Hexahedron.hxx
   StdMeshers_Cartesian_3D.hxx
   StdMeshers_Cartesian_VL.hxx
   StdMeshers_QuadFromMedialAxis_1D2D.hxx
@@ -184,6 +186,8 @@ SET(StdMeshers_SOURCES
   StdMeshers_ViscousLayers2D.cxx
   StdMeshers_Projection_1D2D.cxx
   StdMeshers_CartesianParameters3D.cxx
+  StdMeshers_Cartesian_3D_Grid.cxx
+  StdMeshers_Cartesian_3D_Hexahedron.cxx
   StdMeshers_Cartesian_3D.cxx
   StdMeshers_Cartesian_VL.cxx
   StdMeshers_Adaptive1D.cxx

src/StdMeshers/StdMeshers_CartesianParameters3D.cxx

@@ -361,9 +361,13 @@ void StdMeshers_CartesianParameters3D::ComputeCoordinates(const double    x0,
         ++iCell;
       }
     }
-    const double lastCellLen = coords.back() - coords[ coords.size() - 2 ];
-    if ( fabs( coords.back() - p1 ) > 0.5 * lastCellLen )
+    if (coords.size() < 2)
       coords.push_back ( p1 );
+    else {
+      const double lastCellLen = coords.back() - coords[ coords.size() - 2 ];
+      if ( fabs( coords.back() - p1 ) > 0.5 * lastCellLen )
+        coords.push_back ( p1 );
+    }
   }
 
   // correct coords if a forced point is too close to a neighbor node

src/StdMeshers/StdMeshers_Cartesian_3D_Grid.hxx

@@ -0,0 +1,702 @@
+// Copyright (C) 2016-2024  CEA, EDF
+//
+// This library is free software; you can redistribute it and/or
+// modify it under the terms of the GNU Lesser General Public
+// License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+//
+// This library is distributed in the hope that it will be useful,
+// but WITHOUT ANY WARRANTY; without even the implied warranty of
+// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+// Lesser General Public License for more details.
+//
+// You should have received a copy of the GNU Lesser General Public
+// License along with this library; if not, write to the Free Software
+// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
+//
+// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
+//
+//  File   : StdMeshers_Cartesian_3D_Grid.hxx
+//  Module : SMESH
+//  Purpose: Make BodyFitting mesh algorithm more modular and testable
+//
+
+#ifndef _SMESH_Cartesian_3D_GRID_HXX_
+#define _SMESH_Cartesian_3D_GRID_HXX_
+
+#ifdef WITH_TBB
+#include <tbb/parallel_for.h>
+#endif
+
+#include <utilities.h>
+
+// STD
+#include <algorithm>
+#include <map>
+#include <vector>
+#include <memory>
+#include <mutex>
+#include <thread>
+
+// SMESH
+#include "SMESH_StdMeshers.hxx"
+#include "StdMeshers_FaceSide.hxx"
+#include "StdMeshers_CartesianParameters3D.hxx"
+
+#include <ObjectPool.hxx>
+#include <SMDS_LinearEdge.hxx>
+#include <SMDS_MeshNode.hxx>
+#include <SMDS_VolumeOfNodes.hxx>
+#include <SMDS_VolumeTool.hxx>
+#include <SMESHDS_Mesh.hxx>
+#include <SMESH_Block.hxx>
+#include <SMESH_Comment.hxx>
+#include <SMESH_ControlsDef.hxx>
+#include <SMESH_Mesh.hxx>
+#include <SMESH_MeshAlgos.hxx>
+#include <SMESH_MeshEditor.hxx>
+#include <SMESH_MesherHelper.hxx>
+#include <SMESH_subMesh.hxx>
+#include <SMESH_subMeshEventListener.hxx>
+
+#include <utilities.h>
+#include <Utils_ExceptHandlers.hxx>
+
+#include <GEOMUtils.hxx>
+
+//OCC
+#include <BRepAdaptor_Curve.hxx>
+#include <BRepAdaptor_Surface.hxx>
+#include <BRepBndLib.hxx>
+#include <BRepBuilderAPI_Copy.hxx>
+#include <BRepBuilderAPI_MakeFace.hxx>
+#include <BRepTools.hxx>
+#include <BRepTopAdaptor_FClass2d.hxx>
+#include <BRep_Builder.hxx>
+#include <BRep_Tool.hxx>
+#include <Bnd_B3d.hxx>
+#include <Bnd_Box.hxx>
+#include <ElSLib.hxx>
+#include <GCPnts_UniformDeflection.hxx>
+#include <Geom2d_BSplineCurve.hxx>
+#include <Geom2d_BezierCurve.hxx>
+#include <Geom2d_TrimmedCurve.hxx>
+#include <GeomAPI_ProjectPointOnSurf.hxx>
+#include <GeomLib.hxx>
+#include <Geom_BSplineCurve.hxx>
+#include <Geom_BSplineSurface.hxx>
+#include <Geom_BezierCurve.hxx>
+#include <Geom_BezierSurface.hxx>
+#include <Geom_RectangularTrimmedSurface.hxx>
+#include <Geom_TrimmedCurve.hxx>
+#include <IntAna_IntConicQuad.hxx>
+#include <IntAna_IntLinTorus.hxx>
+#include <IntAna_Quadric.hxx>
+#include <IntCurveSurface_TransitionOnCurve.hxx>
+#include <IntCurvesFace_Intersector.hxx>
+#include <Poly_Triangulation.hxx>
+#include <Precision.hxx>
+#include <TopExp.hxx>
+#include <TopExp_Explorer.hxx>
+#include <TopLoc_Location.hxx>
+#include <TopTools_DataMapOfShapeInteger.hxx>
+#include <TopTools_IndexedMapOfShape.hxx>
+#include <TopTools_MapOfShape.hxx>
+#include <TopoDS.hxx>
+#include <TopoDS_Compound.hxx>
+#include <TopoDS_Face.hxx>
+#include <TopoDS_TShape.hxx>
+#include <gp_Cone.hxx>
+#include <gp_Cylinder.hxx>
+#include <gp_Lin.hxx>
+#include <gp_Pln.hxx>
+#include <gp_Pnt2d.hxx>
+#include <gp_Sphere.hxx>
+#include <gp_Torus.hxx>
+
+// All utility structs used in Grid and hexahedron class will be included here
+// Ideally each one of this should define their own testable class
+namespace StdMeshers
+{
+namespace Cartesian3D
+{
+  typedef int                     TGeomID; // IDs of sub-shapes
+  typedef TopTools_ShapeMapHasher TShapeHasher; // non-oriented shape hasher
+  typedef std::array< int, 3 >    TIJK;
+
+  typedef std::map< TGeomID, std::vector< TGeomID > > TEdge2faceIDsMap;
+
+  const TGeomID theUndefID = 1e+9;
+
+  //=============================================================================
+  // Definitions of internal utils
+  // --------------------------------------------------------------------------
+  enum Transition {
+    Trans_TANGENT = IntCurveSurface_Tangent,
+    Trans_IN      = IntCurveSurface_In,
+    Trans_OUT     = IntCurveSurface_Out,
+    Trans_APEX,
+    Trans_INTERNAL // for INTERNAL FACE
+  };
+
+  // --------------------------------------------------------------------------
+  /*!
+   * \brief Sub-entities of a FACE neighboring its concave VERTEX.
+   *        Help to avoid linking nodes on EDGEs that seem connected
+   *        by the concave FACE but the link actually lies outside the FACE
+   */
+  struct ConcaveFace
+  {
+    TGeomID _concaveFace;
+    TGeomID _edge1, _edge2;
+    TGeomID _v1,    _v2;
+    ConcaveFace( int f=0, int e1=0, int e2=0, int v1=0, int v2=0 )
+      : _concaveFace(f), _edge1(e1), _edge2(e2), _v1(v1), _v2(v2) {}
+    bool HasEdge( TGeomID edge ) const { return edge == _edge1 || edge == _edge2; }
+    bool HasVertex( TGeomID v  ) const { return v == _v1 || v == _v2; }
+    void SetEdge( TGeomID edge ) { ( _edge1 ? _edge2 : _edge1 ) = edge; }
+    void SetVertex( TGeomID v  ) { ( _v1 ? _v2 : _v1 ) = v; }
+  };
+  typedef NCollection_DataMap< TGeomID, ConcaveFace > TConcaveVertex2Face;
+
+  // --------------------------------------------------------------------------
+  /*!
+   * \brief Container of IDs of SOLID sub-shapes
+   */
+  class Solid // sole SOLID contains all sub-shapes
+  {
+    TGeomID             _id; // SOLID id
+    bool                _hasInternalFaces;
+    TConcaveVertex2Face _concaveVertex; // concave VERTEX -> ConcaveFace
+  public:
+    virtual ~Solid() {}
+    virtual bool Contains( TGeomID /*subID*/ ) const { return true; }
+    virtual bool ContainsAny( const std::vector< TGeomID>& /*subIDs*/ ) const { return true; }
+    virtual TopAbs_Orientation Orientation( const TopoDS_Shape& s ) const { return s.Orientation(); }
+    virtual bool IsOutsideOriented( TGeomID /*faceID*/ ) const { return true; }
+    void SetID( TGeomID id ) { _id = id; }
+    TGeomID ID() const { return _id; }
+    void SetHasInternalFaces( bool has ) { _hasInternalFaces = has; }
+    bool HasInternalFaces() const { return _hasInternalFaces; }
+    void SetConcave( TGeomID V, TGeomID F, TGeomID E1, TGeomID E2, TGeomID V1, TGeomID V2  )
+    { _concaveVertex.Bind( V, ConcaveFace{ F, E1, E2, V1, V2 }); }
+    bool HasConcaveVertex() const { return !_concaveVertex.IsEmpty(); }
+    const ConcaveFace* GetConcave( TGeomID V ) const { return _concaveVertex.Seek( V ); }
+  };
+
+  // --------------------------------------------------------------------------
+  class OneOfSolids : public Solid
+  {
+    TColStd_MapOfInteger _subIDs;
+    TopTools_MapOfShape  _faces; // keep FACE orientation
+    TColStd_MapOfInteger _outFaceIDs; // FACEs of shape_to_mesh oriented outside the SOLID
+  public:
+    void Init( const TopoDS_Shape& solid,
+               TopAbs_ShapeEnum    subType,
+               const SMESHDS_Mesh* mesh );
+    virtual bool Contains( TGeomID i ) const { return i == ID() || _subIDs.Contains( i ); }
+    virtual bool ContainsAny( const std::vector< TGeomID>& subIDs ) const
+    {
+      for ( size_t i = 0; i < subIDs.size(); ++i ) if ( Contains( subIDs[ i ])) return true;
+      return false;
+    }
+    virtual TopAbs_Orientation Orientation( const TopoDS_Shape& face ) const
+    {
+      const TopoDS_Shape& sInMap = const_cast< OneOfSolids* >(this)->_faces.Added( face );
+      return sInMap.Orientation();
+    }
+    virtual bool IsOutsideOriented( TGeomID faceID ) const
+    {
+      return faceID == 0 || _outFaceIDs.Contains( faceID );
+    }
+  };
+
+  // --------------------------------------------------------------------------
+  /*!
+   * \brief Hold a vector of TGeomID and clear it at destruction
+   */
+  class GeomIDVecHelder
+  {
+    typedef std::vector< TGeomID > TVector;
+    const TVector& myVec;
+    bool           myOwn;
+
+  public:
+    GeomIDVecHelder( const TVector& idVec, bool isOwner ): myVec( idVec ), myOwn( isOwner ) {}
+    GeomIDVecHelder( const GeomIDVecHelder& holder ): myVec( holder.myVec ), myOwn( holder.myOwn )
+    {
+      const_cast< bool& >( holder.myOwn ) = false;
+    }
+    ~GeomIDVecHelder() { if ( myOwn ) const_cast<TVector&>( myVec ).clear(); }
+    size_t size() const { return myVec.size(); }
+    TGeomID operator[]( size_t i ) const { return i < size() ? myVec[i] : theUndefID; }
+    bool operator==( const GeomIDVecHelder& other ) const { return myVec == other.myVec; }
+    bool contain( const TGeomID& id ) const {
+      return std::find( myVec.begin(), myVec.end(), id ) != myVec.end();
+    }
+    TGeomID otherThan( const TGeomID& id ) const {
+      for ( const TGeomID& id2 : myVec )
+        if ( id != id2 )
+          return id2;
+      return theUndefID;
+    }
+    TGeomID oneCommon( const GeomIDVecHelder& other ) const {
+      TGeomID common = theUndefID;
+      for ( const TGeomID& id : myVec )
+        if ( other.contain( id ))
+        {
+          if ( common != theUndefID )
+            return theUndefID;
+          common = id;
+        }
+      return common;
+    }
+  };
+
+  // --------------------------------------------------------------------------
+  /*!
+   * \brief Geom data
+   */
+  struct Geometry
+  {
+    TopoDS_Shape                _mainShape;
+    std::vector< std::vector< TGeomID > > _solidIDsByShapeID;// V/E/F ID -> SOLID IDs
+    Solid                            _soleSolid;
+    std::map< TGeomID, OneOfSolids > _solidByID;
+    TColStd_MapOfInteger        _boundaryFaces; // FACEs on boundary of mesh->ShapeToMesh()
+    TColStd_MapOfInteger        _strangeEdges; // EDGEs shared by strange FACEs
+    TGeomID                     _extIntFaceID; // pseudo FACE - extension of INTERNAL FACE
+
+    TopTools_DataMapOfShapeInteger _shape2NbNodes; // nb of pre-existing nodes on shapes
+
+    SMESH::Controls::ElementsOnShape _edgeClassifier;
+    SMESH::Controls::ElementsOnShape _vertexClassifier;
+
+    bool IsOneSolid() const { return _solidByID.size() < 2; }
+    GeomIDVecHelder GetSolidIDsByShapeID( const std::vector< TGeomID >& shapeIDs ) const;
+  };
+
+  // --------------------------------------------------------------------------
+  /*!
+   * \brief Common data of any intersection between a Grid and a shape
+   */
+  struct B_IntersectPoint
+  {
+    // This two class members are being updated in a non thread safe way.
+    // See Add method modify _node and _faceIDs class members dinamicaly during execution
+    // of Hexahedron.compute() method.
+    // std::mutex _mutex;
+    mutable const SMDS_MeshNode*   _node;
+    mutable std::vector< TGeomID > _faceIDs;
+
+    B_IntersectPoint(): _node(NULL) {}
+    bool Add( const std::vector< TGeomID >& fIDs, const SMDS_MeshNode* n=NULL ) const;
+    TGeomID HasCommonFace( const B_IntersectPoint * other, TGeomID avoidFace=-1 ) const;
+    size_t GetCommonFaces( const B_IntersectPoint * other, TGeomID * commonFaces ) const;
+    bool IsOnFace( TGeomID faceID ) const;
+    virtual ~B_IntersectPoint() {}
+  };
+
+  // --------------------------------------------------------------------------
+  /*!
+   * \brief Data of intersection between a GridLine and a TopoDS_Face
+   */
+  struct F_IntersectPoint : public B_IntersectPoint
+  {
+    double             _paramOnLine;
+    double             _u, _v;
+    mutable Transition _transition;
+    mutable size_t     _indexOnLine;
+
+    bool operator< ( const F_IntersectPoint& o ) const { 
+      return _paramOnLine < o._paramOnLine;         
+    }
+  };
+
+  // --------------------------------------------------------------------------
+  /*!
+   * \brief Data of intersection between GridPlanes and a TopoDS_EDGE
+   */
+  struct E_IntersectPoint : public B_IntersectPoint
+  {
+    gp_Pnt  _point;
+    double  _uvw[3];
+    TGeomID _shapeID; // ID of EDGE or VERTEX
+  };
+
+  // --------------------------------------------------------------------------
+  /*!
+    * \brief A line of the grid and its intersections with 2D geometry
+    */
+  struct GridLine
+  {
+    gp_Lin _line;
+    double _length; // line length
+    std::multiset< F_IntersectPoint > _intPoints;
+
+    void RemoveExcessIntPoints( const double tol );
+    TGeomID GetSolidIDBefore( std::multiset< F_IntersectPoint >::iterator ip,
+                              const TGeomID                          prevID,
+                              const Geometry&                        geom);
+  };
+  // --------------------------------------------------------------------------
+  /*!
+   * \brief Planes of the grid used to find intersections of an EDGE with a hexahedron
+   */
+  struct GridPlanes
+  {
+    gp_XYZ                _zNorm;
+    std::vector< gp_XYZ > _origins; // origin points of all planes in one direction
+    std::vector< double > _zProjs;  // projections of origins to _zNorm
+  };
+  // --------------------------------------------------------------------------
+  /*!
+   * \brief Iterator on the parallel grid lines of one direction
+   */
+  struct LineIndexer
+  {
+    size_t _size  [3];
+    size_t _curInd[3];
+    size_t _iVar1, _iVar2, _iConst;
+    std::string _name1, _name2, _nameConst;
+    LineIndexer() {}
+    LineIndexer( size_t sz1, size_t sz2, size_t sz3,
+                 size_t iv1, size_t iv2, size_t iConst,
+                 const std::string& nv1, const std::string& nv2, const std::string& nConst )
+    {
+      _size[0] = sz1; _size[1] = sz2; _size[2] = sz3;
+      _curInd[0] = _curInd[1] = _curInd[2] = 0;
+      _iVar1 = iv1; _iVar2 = iv2; _iConst = iConst;
+      _name1 = nv1; _name2 = nv2; _nameConst = nConst;
+    }
+
+    size_t I() const { return _curInd[0]; }
+    size_t J() const { return _curInd[1]; }
+    size_t K() const { return _curInd[2]; }
+    void SetIJK( size_t i, size_t j, size_t k )
+    {
+      _curInd[0] = i; _curInd[1] = j; _curInd[2] = k;
+    }
+    void SetLineIndex(size_t i)
+    {
+      _curInd[_iVar2] = i / _size[_iVar1];
+      _curInd[_iVar1] = i % _size[_iVar1];
+    }
+    void operator++()
+    {
+      if ( ++_curInd[_iVar1] == _size[_iVar1] )
+        _curInd[_iVar1] = 0, ++_curInd[_iVar2];
+    }
+    bool More() const { return _curInd[_iVar2] < _size[_iVar2]; }
+    size_t LineIndex   () const { return _curInd[_iVar1] + _curInd[_iVar2]* _size[_iVar1]; }
+    size_t LineIndex10 () const { return (_curInd[_iVar1] + 1 ) + _curInd[_iVar2]* _size[_iVar1]; }
+    size_t LineIndex01 () const { return _curInd[_iVar1] + (_curInd[_iVar2] + 1 )* _size[_iVar1]; }
+    size_t LineIndex11 () const { return (_curInd[_iVar1] + 1 ) + (_curInd[_iVar2] + 1 )* _size[_iVar1]; }
+    void SetIndexOnLine (size_t i)  { _curInd[ _iConst ] = i; }
+    bool IsValidIndexOnLine (size_t i) const { return  i < _size[ _iConst ]; }
+    size_t NbLines() const { return _size[_iVar1] * _size[_iVar2]; }
+  };
+
+  class Tools 
+  { 
+    public:
+      Tools() = delete;
+      
+      //================================================================================
+      /*!
+      * \brief computes exact bounding box with axes parallel to given ones
+      */
+      //================================================================================
+      static void GetExactBndBox( const std::vector< TopoDS_Shape >& faceVec, const double* axesDirs, Bnd_Box& shapeBox );
+  };
+
+  class STDMESHERS_EXPORT Grid
+  {
+    public:
+    std::vector< double >   _coords[3]; // coordinates of grid nodes
+    gp_XYZ                  _axes  [3]; // axis directions
+    std::vector< GridLine > _lines [3]; //    in 3 directions
+    double                  _tol, _minCellSize;
+    gp_XYZ                 _origin;
+    gp_Mat                 _invB; // inverted basis of _axes
+
+    // index shift within _nodes of nodes of a cell from the 1st node
+    int                    _nodeShift[8];
+
+    std::vector< const SMDS_MeshNode* >    _nodes;          // mesh nodes at grid nodes
+    std::vector< const SMDS_MeshNode* >    _allBorderNodes; // mesh nodes between the bounding box and the geometry boundary
+
+    std::vector< const F_IntersectPoint* > _gridIntP; // grid node intersection with geometry
+    ObjectPool< E_IntersectPoint >        _edgeIntPool; // intersections with EDGEs
+    ObjectPool< F_IntersectPoint >        _extIntPool; // intersections with extended INTERNAL FACEs
+    //list< E_IntersectPoint >          _edgeIntP; // intersections with EDGEs
+
+    Geometry                          _geometry;
+    bool                              _toAddEdges;
+    bool                              _toCreateFaces;
+    bool                              _toConsiderInternalFaces;
+    bool                              _toUseThresholdForInternalFaces;
+    double                            _sizeThreshold;
+    bool                              _toUseQuanta;
+    double                            _quanta;
+
+    SMESH_MesherHelper*               _helper;
+
+    size_t CellIndex( size_t i, size_t j, size_t k ) const
+    {
+      return i + j*(_coords[0].size()-1) + k*(_coords[0].size()-1)*(_coords[1].size()-1);
+    }
+    size_t NodeIndex( size_t i, size_t j, size_t k ) const
+    {
+      return i + j*_coords[0].size() + k*_coords[0].size()*_coords[1].size();
+    }
+    size_t NodeIndex( const TIJK& ijk ) const
+    {
+      return NodeIndex( ijk[0], ijk[1], ijk[2] );
+    }
+    size_t NodeIndexDX() const { return 1; }
+    size_t NodeIndexDY() const { return _coords[0].size(); }
+    size_t NodeIndexDZ() const { return _coords[0].size() * _coords[1].size(); }
+
+    LineIndexer GetLineIndexer(size_t iDir) const;
+    size_t GetLineDir( const GridLine* line, size_t & index ) const;
+
+    E_IntersectPoint* Add( const E_IntersectPoint& ip )
+    {
+      E_IntersectPoint* eip = _edgeIntPool.getNew();
+      *eip = ip;
+      return eip;
+    }
+    void Remove( E_IntersectPoint* eip ) { _edgeIntPool.destroy( eip ); }
+
+    TGeomID ShapeID( const TopoDS_Shape& s ) const;
+    const TopoDS_Shape& Shape( TGeomID id ) const;
+    TopAbs_ShapeEnum ShapeType( TGeomID id ) const { return Shape(id).ShapeType(); }
+    void InitGeometry( const TopoDS_Shape& theShape );
+    void InitClassifier( const TopoDS_Shape&        mainShape,
+                        TopAbs_ShapeEnum           shapeType,
+                        SMESH::Controls::ElementsOnShape& classifier );
+    void GetEdgesToImplement( std::map< TGeomID, std::vector< TGeomID > > & edge2faceMap,
+                              const TopoDS_Shape&                 shape,
+                              const std::vector< TopoDS_Shape >&       faces );
+    void SetSolidFather( const TopoDS_Shape& s, const TopoDS_Shape& theShapeToMesh );
+    bool IsShared( TGeomID faceID ) const;
+    bool IsAnyShared( const std::vector< TGeomID >& faceIDs ) const;
+    bool IsInternal( TGeomID faceID ) const {
+      return ( faceID == PseudoIntExtFaceID() ||
+              Shape( faceID ).Orientation() == TopAbs_INTERNAL ); }
+    bool IsSolid( TGeomID shapeID ) const {
+      if ( _geometry.IsOneSolid() ) return _geometry._soleSolid.ID() == shapeID;
+      else                          return _geometry._solidByID.count( shapeID ); }
+    bool IsStrangeEdge( TGeomID id ) const { return _geometry._strangeEdges.Contains( id ); }
+    TGeomID PseudoIntExtFaceID() const { return _geometry._extIntFaceID; }
+    Solid* GetSolid( TGeomID solidID = 0 );
+    Solid* GetOneOfSolids( TGeomID solidID );
+    const std::vector< TGeomID > & GetSolidIDs( TGeomID subShapeID ) const;
+    bool IsCorrectTransition( TGeomID faceID, const Solid* solid );
+    bool IsBoundaryFace( TGeomID face ) const { return _geometry._boundaryFaces.Contains( face ); }
+    void SetOnShape( const SMDS_MeshNode* n, const F_IntersectPoint& ip,
+                    TopoDS_Vertex* vertex = nullptr, bool unset = false );
+    void UpdateFacesOfVertex( const B_IntersectPoint& ip, const TopoDS_Vertex& vertex );
+    bool IsToCheckNodePos() const { return !_toAddEdges && _toCreateFaces; }
+    bool IsToRemoveExcessEntities() const { return !_toAddEdges; }
+
+    void SetCoordinates(const std::vector<double>& xCoords,
+                        const std::vector<double>& yCoords,
+                        const std::vector<double>& zCoords,
+                        const double*         axesDirs,
+                        const Bnd_Box&        bndBox );
+    void ComputeUVW(const gp_XYZ& p, double uvw[3]);
+    void ComputeNodes(SMESH_MesherHelper& helper);
+    bool GridInitAndInterserctWithShape( const TopoDS_Shape& theShape,
+                                         std::map< TGeomID, std::vector< TGeomID > >& edge2faceIDsMap, 
+                                         const StdMeshers_CartesianParameters3D* hyp,
+                                         const int numOfThreads,
+                                         bool computeCanceled );
+  };
+
+  // Implement parallel computation of Hexa with c++ thread implementation
+  template<typename Iterator, class Function>
+  void parallel_for(const Iterator& first, const Iterator& last, Function&& f, const unsigned int nthreads = 1)
+  {
+    MESSAGE("Start parallel computation of Hexa with c++ threads...");
+
+    assert(nthreads > 0);
+
+    const unsigned int numTasksTotal = last - first;
+    std::vector<std::thread> threads;
+    Iterator it = first;
+
+    MESSAGE("Number of elements to compute: " << numTasksTotal << "; num of threads: " << nthreads);
+
+    // Distribute tasks among threads
+    if (numTasksTotal <= nthreads)
+    {
+      // A simple case - just one task executed in one thread.
+      // TODO: check if it's faster to do it sequentially
+      threads.reserve(numTasksTotal);
+      for (; it < last; ++it)
+      {
+        threads.emplace_back(f, std::ref(*it));
+      }
+
+      // This line for debug in sequential mode
+      // std::for_each(it, last, f);
+    }
+    else
+    {
+      // Calculate how to distribute elements among threads evenly
+      const unsigned int numTasksInThread = numTasksTotal / nthreads;
+      MESSAGE("Number of tasks in thread: " << numTasksInThread);
+
+      // Store the numbers of tasks per thread
+      std::vector<unsigned int> distTasksInThreads(nthreads, numTasksInThread);
+
+      // Distribute a remainder among saved numbers
+      const unsigned int remainder = numTasksTotal % nthreads;
+      MESSAGE("Remainder of tasks " << remainder << " will be evenly distributed among threads");
+      for (unsigned int i = 0; i < remainder; ++i)
+      {
+        ++distTasksInThreads[i];
+      }
+
+      // Create threads for each number of tasks
+      threads.reserve(nthreads);
+      for (const auto i : distTasksInThreads)
+      {
+        Iterator curLast = it + i;
+
+        // Pass iterators by value and the function by reference!
+        auto lambda = [=,&f](){ std::for_each(it, curLast, f); };
+
+        // Create a thread
+        threads.emplace_back(lambda);
+
+        // Advance iterator to the next step
+        it = curLast;
+      }
+    }
+
+    std::for_each(threads.begin(), threads.end(), [](std::thread& x){ x.join(); });
+
+    MESSAGE("Parallel computation was finished successfully");
+  }
+
+  // --------------------------------------------------------------------------
+  /*!
+   * \brief Intersector of TopoDS_Face with all GridLine's
+   */
+  struct FaceGridIntersector
+  {
+    TopoDS_Face _face;
+    TGeomID     _faceID;
+    Grid*       _grid;
+    Bnd_Box     _bndBox;
+    IntCurvesFace_Intersector* _surfaceInt;
+    std::vector< std::pair< GridLine*, F_IntersectPoint > > _intersections;
+
+    FaceGridIntersector(): _grid(0), _surfaceInt(0) {}
+    void Intersect();
+
+    void StoreIntersections()
+    {
+      for ( size_t i = 0; i < _intersections.size(); ++i )
+      {
+        std::multiset< F_IntersectPoint >::iterator ip =
+          _intersections[i].first->_intPoints.insert( _intersections[i].second );
+        ip->_faceIDs.reserve( 1 );
+        ip->_faceIDs.push_back( _faceID );
+      }
+    }
+    const Bnd_Box& GetFaceBndBox()
+    {
+      GetCurveFaceIntersector();
+      return _bndBox;
+    }
+    IntCurvesFace_Intersector* GetCurveFaceIntersector()
+    {
+      if ( !_surfaceInt )
+      {
+        _surfaceInt = new IntCurvesFace_Intersector( _face, Precision::PConfusion() );
+        _bndBox     = _surfaceInt->Bounding();
+        if ( _bndBox.IsVoid() )
+          BRepBndLib::Add (_face, _bndBox);
+      }
+      return _surfaceInt;
+    }
+#ifdef WITH_TBB    
+    bool IsThreadSafe(std::set< const Standard_Transient* >& noSafeTShapes) const;
+#endif
+  };
+  
+#ifdef WITH_TBB
+  // --------------------------------------------------------------------------
+  /*!
+   * \brief Structure intersecting certain nb of faces with GridLine's in one thread
+   */
+  struct ParallelIntersector
+  {
+    std::vector< FaceGridIntersector >& _faceVec;
+    ParallelIntersector( std::vector< FaceGridIntersector >& faceVec): _faceVec(faceVec){}
+    void operator() ( const tbb::blocked_range<size_t>& r ) const
+    {
+      for ( size_t i = r.begin(); i != r.end(); ++i )
+        _faceVec[i].Intersect();
+    }
+  };
+#endif
+
+  template<typename Type>
+  void computeGridIntersection( Type faceGridIntersector )
+  {
+    faceGridIntersector.Intersect();
+  }
+
+  // --------------------------------------------------------------------------
+  /*!
+   * \brief Intersector of a surface with a GridLine
+   */
+  struct FaceLineIntersector
+  {
+    double      _tol;
+    double      _u, _v, _w; // params on the face and the line
+    Transition  _transition; // transition at intersection (see IntCurveSurface.cdl)
+    Transition  _transIn, _transOut; // IN and OUT transitions depending of face orientation
+
+    gp_Pln      _plane;
+    gp_Cylinder _cylinder;
+    gp_Cone     _cone;
+    gp_Sphere   _sphere;
+    gp_Torus    _torus;
+    IntCurvesFace_Intersector* _surfaceInt;
+
+    std::vector< F_IntersectPoint > _intPoints;
+
+    void IntersectWithPlane   (const GridLine& gridLine);
+    void IntersectWithCylinder(const GridLine& gridLine);
+    void IntersectWithCone    (const GridLine& gridLine);
+    void IntersectWithSphere  (const GridLine& gridLine);
+    void IntersectWithTorus   (const GridLine& gridLine);
+    void IntersectWithSurface (const GridLine& gridLine);
+
+    /*
+     * Return true if (_u,_v) is on the face
+     */
+    bool UVIsOnFace() const
+    {
+      TopAbs_State state = _surfaceInt->ClassifyUVPoint(gp_Pnt2d( _u,_v ));
+      return ( state == TopAbs_IN || state == TopAbs_ON );
+    }
+    void addIntPoint(const bool toClassify=true);
+    bool isParamOnLineOK( const double linLength )
+    {
+      return -_tol < _w && _w < linLength + _tol;
+    }
+    FaceLineIntersector():_surfaceInt(0) {}
+    ~FaceLineIntersector() { if (_surfaceInt ) delete _surfaceInt; _surfaceInt = 0; }
+  };
+
+} // end namespace Cartesian3D
+} // end namespace StdMeshers
+
+#endif

src/StdMeshers/StdMeshers_Cartesian_3D_Hexahedron.hxx

@@ -0,0 +1,653 @@
+// Copyright (C) 2016-2024  CEA, EDF
+//
+// This library is free software; you can redistribute it and/or
+// modify it under the terms of the GNU Lesser General Public
+// License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+//
+// This library is distributed in the hope that it will be useful,
+// but WITHOUT ANY WARRANTY; without even the implied warranty of
+// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+// Lesser General Public License for more details.
+//
+// You should have received a copy of the GNU Lesser General Public
+// License along with this library; if not, write to the Free Software
+// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
+//
+// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
+//
+//  File   : StdMeshers_Cartesian_3D_Hexahedron.hxx
+//  Module : SMESH
+//  Purpose: Make BodyFitting mesh algorithm more modular and testable
+//
+
+#ifndef _SMESH_Cartesian_3D_HEXAHEDRON_HXX_
+#define _SMESH_Cartesian_3D_HEXAHEDRON_HXX_
+
+// BOOST
+#include <boost/container/flat_map.hpp>
+
+// STD
+#include <utilities.h>
+#include <vector>
+
+// SMESH
+#include "SMESH_StdMeshers.hxx"
+#include "StdMeshers_Cartesian_3D_Grid.hxx"
+
+namespace StdMeshers
+{
+namespace Cartesian3D
+{
+  // --------------------------------------------------------------------------
+  /*!
+   * \brief Return cells sharing a link
+   */
+  struct CellsAroundLink
+  {
+    int    _iDir;
+    int    _dInd[4][3];
+    size_t _nbCells[3];
+    int    _i,_j,_k;
+    StdMeshers::Cartesian3D::Grid*  _grid;
+
+    CellsAroundLink( StdMeshers::Cartesian3D::Grid* grid, int iDir ):
+      _iDir( iDir ),
+      _dInd{ {0,0,0}, {0,0,0}, {0,0,0}, {0,0,0} },
+      _nbCells{ grid->_coords[0].size() - 1,
+      grid->_coords[1].size() - 1,
+      grid->_coords[2].size() - 1 },
+      _grid( grid )
+    {
+      const int iDirOther[3][2] = {{ 1,2 },{ 0,2 },{ 0,1 }};
+      _dInd[1][ iDirOther[iDir][0] ] = -1;
+      _dInd[2][ iDirOther[iDir][1] ] = -1;
+      _dInd[3][ iDirOther[iDir][0] ] = -1; _dInd[3][ iDirOther[iDir][1] ] = -1;
+    }
+    void Init( int i, int j, int k, int link12 = 0 )
+    {
+      int iL = link12 % 4;
+      _i = i - _dInd[iL][0];
+      _j = j - _dInd[iL][1];
+      _k = k - _dInd[iL][2];
+    }
+    bool GetCell( int iL, int& i, int& j, int& k, int& cellIndex, int& linkIndex )
+    {
+      i =  _i + _dInd[iL][0];
+      j =  _j + _dInd[iL][1];
+      k =  _k + _dInd[iL][2];
+      if ( i < 0 || i >= (int)_nbCells[0] ||
+           j < 0 || j >= (int)_nbCells[1] ||
+           k < 0 || k >= (int)_nbCells[2] )
+        return false;
+      cellIndex = _grid->CellIndex( i,j,k );
+      linkIndex = iL + _iDir * 4;
+      return true;
+    }
+  };
+
+  // --------------------------------------------------------------------------
+  /*!
+   * \brief Class representing topology of the hexahedron and creating a mesh
+   *        volume basing on analysis of hexahedron intersection with geometry
+   */
+  class STDMESHERS_EXPORT Hexahedron
+  {
+    // --------------------------------------------------------------------------------
+    struct _Face;
+    struct _Link;
+    enum IsInternalFlag { IS_NOT_INTERNAL, IS_INTERNAL, IS_CUT_BY_INTERNAL_FACE };
+    // --------------------------------------------------------------------------------
+    struct _Node //!< node either at a hexahedron corner or at intersection
+    {
+      const SMDS_MeshNode*    _node;        // mesh node at hexahedron corner
+      const SMDS_MeshNode*    _boundaryCornerNode; // missing mesh node due to hex truncation on the boundary
+      const StdMeshers::Cartesian3D::B_IntersectPoint* _intPoint;
+      const _Face*            _usedInFace;
+      char                    _isInternalFlags;
+
+      _Node(const SMDS_MeshNode* n=0, const StdMeshers::Cartesian3D::B_IntersectPoint* ip=0)
+        :_node(n), _boundaryCornerNode(0), _intPoint(ip), _usedInFace(0), _isInternalFlags(0) {}
+      const SMDS_MeshNode*    Node() const
+      { return ( _intPoint && _intPoint->_node ) ? _intPoint->_node : _node; }
+      const SMDS_MeshNode*    BoundaryNode() const
+      { return _node ? _node : _boundaryCornerNode; }
+      const StdMeshers::Cartesian3D::E_IntersectPoint* EdgeIntPnt() const
+      { return static_cast< const StdMeshers::Cartesian3D::E_IntersectPoint* >( _intPoint ); }
+      const StdMeshers::Cartesian3D::F_IntersectPoint* FaceIntPnt() const
+      { return static_cast< const StdMeshers::Cartesian3D::F_IntersectPoint* >( _intPoint ); }
+      const std::vector< StdMeshers::Cartesian3D::TGeomID >& faces() const { return _intPoint->_faceIDs; }
+      StdMeshers::Cartesian3D::TGeomID face(size_t i) const { return _intPoint->_faceIDs[ i ]; }
+      void SetInternal( IsInternalFlag intFlag ) { _isInternalFlags |= intFlag; }
+      bool IsCutByInternal() const { return _isInternalFlags & IS_CUT_BY_INTERNAL_FACE; }
+      bool IsUsedInFace( const _Face* polygon = 0 )
+      {
+        return polygon ? ( _usedInFace == polygon ) : bool( _usedInFace );
+      }
+      StdMeshers::Cartesian3D::TGeomID IsLinked( const StdMeshers::Cartesian3D::B_IntersectPoint* other,
+                                                 StdMeshers::Cartesian3D::TGeomID avoidFace=-1 ) const // returns id of a common face
+      {
+        return _intPoint ? _intPoint->HasCommonFace( other, avoidFace ) : 0;
+      }
+      bool IsOnFace( StdMeshers::Cartesian3D::TGeomID faceID ) const // returns true if faceID is found
+      {
+        return _intPoint ? _intPoint->IsOnFace( faceID ) : false;
+      }
+      size_t GetCommonFaces( const StdMeshers::Cartesian3D::B_IntersectPoint * other,
+                             StdMeshers::Cartesian3D::TGeomID* common ) const
+      {
+        return _intPoint && other ? _intPoint->GetCommonFaces( other, common ) : 0;
+      }
+      gp_Pnt Point() const
+      {
+        if ( const SMDS_MeshNode* n = Node() )
+          return SMESH_NodeXYZ( n );
+        if ( const StdMeshers::Cartesian3D::E_IntersectPoint* eip =
+             dynamic_cast< const StdMeshers::Cartesian3D::E_IntersectPoint* >( _intPoint ))
+          return eip->_point;
+        return gp_Pnt( 1e100, 0, 0 );
+      }
+      StdMeshers::Cartesian3D::TGeomID ShapeID() const
+      {
+        if ( const StdMeshers::Cartesian3D::E_IntersectPoint* eip =
+             dynamic_cast< const StdMeshers::Cartesian3D::E_IntersectPoint* >( _intPoint ))
+          return eip->_shapeID;
+        return 0;
+      }
+
+      void Add( const StdMeshers::Cartesian3D::E_IntersectPoint* ip );
+      void clear();
+
+      friend std::ostream& operator<<(std::ostream& os, const _Node& node)
+      {
+        if (node._node)
+        {
+          os << "Node at hexahedron corner: ";
+          node._node->Print(os); 
+        }
+        else if (node._intPoint && node._intPoint->_node)
+        {
+          os << "Node at intersection point: ";
+          node._intPoint->_node->Print(os); // intersection point
+        }
+        else
+          os << "mesh node is null\n";
+
+        return os;
+      }
+    };
+
+    // --------------------------------------------------------------------------------
+    struct _Link // link connecting two _Node's
+    {
+      static const std::size_t nodesNum = 2;
+
+      _Node* _nodes[nodesNum];
+      _Face* _faces[nodesNum]; // polygons sharing a link
+      std::vector< const StdMeshers::Cartesian3D::F_IntersectPoint* > _fIntPoints; // GridLine intersections with FACEs
+      std::vector< _Node* > _fIntNodes;   // _Node's at _fIntPoints
+      std::vector< _Link >  _splits;
+      _Link(): _nodes{ 0, 0 }, _faces{ 0, 0 } {}
+
+      void clear();
+
+      friend std::ostream& operator<<(std::ostream& os, const _Link& link)
+      {
+        os << "Link:\n";
+
+        for (std::size_t i = 0; i < nodesNum; ++i)
+        {
+          if (link._nodes[i])
+            os << *link._nodes[i];
+          else
+            os << "link node with index " << i << " is null\n";
+        }
+
+        os << "_fIntPoints: " << link._fIntPoints.size() << '\n';
+        os << "_fIntNodes: " << link._fIntNodes.size() << '\n';
+        os << "_splits: " << link._splits.size() << '\n';
+
+        return os;
+      }
+    };
+
+    // --------------------------------------------------------------------------------
+    struct _OrientedLink
+    {
+      _Link* _link;
+      bool   _reverse;
+      _OrientedLink( _Link* link=0, bool reverse=false ): _link(link), _reverse(reverse) {}
+      void Reverse() { _reverse = !_reverse; }
+      size_t NbResultLinks() const { return _link->_splits.size(); }
+      _OrientedLink ResultLink(int i) const
+      {
+        return _OrientedLink(&_link->_splits[_reverse ? NbResultLinks()-i-1 : i],_reverse);
+      }
+      _Node* FirstNode() const { return _link->_nodes[ _reverse ]; }
+      _Node* LastNode()  const { return _link->_nodes[ !_reverse ]; }
+      operator bool() const { return _link; }
+
+      // returns supporting FACEs
+      std::vector< StdMeshers::Cartesian3D::TGeomID > GetNotUsedFaces
+      (const std::set<StdMeshers::Cartesian3D::TGeomID>& usedIDs ) const
+      {
+        std::vector< StdMeshers::Cartesian3D::TGeomID > faces;
+        const StdMeshers::Cartesian3D::B_IntersectPoint *ip0, *ip1;
+        if (( ip0 = _link->_nodes[0]->_intPoint ) &&
+            ( ip1 = _link->_nodes[1]->_intPoint ))
+        {
+          for ( size_t i = 0; i < ip0->_faceIDs.size(); ++i )
+            if ( ip1->IsOnFace ( ip0->_faceIDs[i] ) &&
+                 !usedIDs.count( ip0->_faceIDs[i] ) )
+              faces.push_back( ip0->_faceIDs[i] );
+        }
+        return faces;
+      }
+
+      bool HasEdgeNodes() const
+      {
+        return ( dynamic_cast< const StdMeshers::Cartesian3D::E_IntersectPoint* >( _link->_nodes[0]->_intPoint ) ||
+                 dynamic_cast< const StdMeshers::Cartesian3D::E_IntersectPoint* >( _link->_nodes[1]->_intPoint ));
+      }
+      int NbFaces() const
+      {
+        return !_link->_faces[0] ? 0 : 1 + bool( _link->_faces[1] );
+      }
+      void AddFace( _Face* f )
+      {
+        if ( _link->_faces[0] ) {
+          _link->_faces[1] = f;
+        }
+        else {
+          _link->_faces[0] = f;
+          _link->_faces[1] = 0;
+        }
+      }
+      void RemoveFace( const _Face* f )
+      {
+        if ( !_link->_faces[0] ) return;
+
+        if ( _link->_faces[1] == f ) {
+          _link->_faces[1] = 0;
+        }
+        else if ( _link->_faces[0] == f ) {
+          _link->_faces[0] = 0;
+          if ( _link->_faces[1] ) {
+            _link->_faces[0] = _link->_faces[1];
+            _link->_faces[1] = 0;
+          }
+        }
+      }
+
+      friend std::ostream& operator<<(std::ostream& os, const _OrientedLink& link)
+      {
+        if (link._link)
+          os << "Oriented " << *link._link;
+        else
+          os << "Oriented link is null\n";
+
+        return os;
+      }
+    };
+
+    // --------------------------------------------------------------------------------
+    struct _SplitIterator //! set to _hexLinks splits on one side of INTERNAL FACEs
+    {
+      struct _Split // data of a link split
+      {
+        int    _linkID;          // hex link ID
+        _Node* _nodes[2];
+        int    _iCheckIteration; // iteration where split is tried as Hexahedron split
+        _Link* _checkedSplit;    // split set to hex links
+        bool   _isUsed;          // used in a volume
+
+        _Split( _Link & split, int iLink ):
+          _linkID( iLink ), _nodes{ split._nodes[0], split._nodes[1] },
+          _iCheckIteration( 0 ), _isUsed( false )
+        {}
+        bool IsCheckedOrUsed( bool used ) const { return used ? _isUsed : _iCheckIteration > 0; }
+      };
+      _Link*                _hexLinks;
+      std::vector< _Split > _splits;
+      int                   _iterationNb;
+      size_t                _nbChecked;
+      size_t                _nbUsed;
+      std::vector< _Node* > _freeNodes; // nodes reached while composing a split set
+
+      _SplitIterator( _Link* hexLinks ):
+        _hexLinks( hexLinks ), _iterationNb(0), _nbChecked(0), _nbUsed(0)
+      {
+        _freeNodes.reserve( 12 );
+        _splits.reserve( 24 );
+        for ( int iL = 0; iL < 12; ++iL )
+          for ( size_t iS = 0; iS < _hexLinks[ iL ]._splits.size(); ++iS )
+            _splits.emplace_back( _hexLinks[ iL ]._splits[ iS ], iL );
+        Next();
+      }
+      bool More() const { return _nbUsed < _splits.size(); }
+      bool Next();
+    };
+
+    // --------------------------------------------------------------------------------
+    struct _Face
+    {
+      SMESH_Block::TShapeID   _name;
+      std::vector< _OrientedLink > _links;       // links on GridLine's
+      std::vector< _Link >         _polyLinks;   // links added to close a polygonal face
+      std::vector< _Node* >        _eIntNodes;   // nodes at intersection with EDGEs
+
+      _Face():_name( SMESH_Block::ID_NONE )
+      {}
+      bool IsPolyLink( const _OrientedLink& ol )
+      {
+        return _polyLinks.empty() ? false :
+          ( &_polyLinks[0] <= ol._link &&  ol._link <= &_polyLinks.back() );
+      }
+      void AddPolyLink(_Node* n0, _Node* n1, _Face* faceToFindEqual=0)
+      {
+        if ( faceToFindEqual && faceToFindEqual != this ) {
+          for ( size_t iL = 0; iL < faceToFindEqual->_polyLinks.size(); ++iL )
+            if ( faceToFindEqual->_polyLinks[iL]._nodes[0] == n1 &&
+                 faceToFindEqual->_polyLinks[iL]._nodes[1] == n0 )
+            {
+              _links.push_back
+                ( _OrientedLink( & faceToFindEqual->_polyLinks[iL], /*reverse=*/true ));
+              return;
+            }
+        }
+        _Link l;
+        l._nodes[0] = n0;
+        l._nodes[1] = n1;
+        _polyLinks.push_back( l );
+        _links.push_back( _OrientedLink( &_polyLinks.back() ));
+      }
+
+      friend std::ostream& operator<<(std::ostream& os, const _Face& face)
+      {
+        os << "Face " << face._name << '\n';
+
+        os << "Links on GridLines: \n";
+        for (const auto& link : face._links)
+        {
+          os << link;
+        }
+
+        os << "Links added to close a polygonal face: \n";
+        for (const auto& link : face._polyLinks)
+        {
+          os << link;
+        }
+
+        os << "Nodes at intersection with EDGEs: \n";
+        for (const auto node : face._eIntNodes)
+        {
+          if (node)
+          {
+            os << *node;
+          }
+        }
+
+        return os;
+      }
+    };
+
+    // --------------------------------------------------------------------------------
+    struct _volumeDef // holder of nodes of a volume mesh element
+    {
+      typedef void* _ptr;
+
+      struct _nodeDef
+      {
+        const SMDS_MeshNode*    _node; // mesh node at hexahedron corner
+        const StdMeshers::Cartesian3D::B_IntersectPoint* _intPoint;
+
+        _nodeDef(): _node(0), _intPoint(0) {}
+        _nodeDef( _Node* n ): _node( n->_node), _intPoint( n->_intPoint ) {}
+        const SMDS_MeshNode*    Node() const
+        { return ( _intPoint && _intPoint->_node ) ? _intPoint->_node : _node; }
+        const StdMeshers::Cartesian3D::E_IntersectPoint* EdgeIntPnt() const
+        { return static_cast< const StdMeshers::Cartesian3D::E_IntersectPoint* >( _intPoint ); }
+        _ptr Ptr() const { return Node() ? (_ptr) Node() : (_ptr) EdgeIntPnt(); }
+        bool operator==(const _nodeDef& other ) const { return Ptr() == other.Ptr(); }
+
+        friend std::ostream& operator<<(std::ostream& os, const _nodeDef& node)
+        {
+          if (node._node)
+          {
+            os << "Node at hexahedron corner: ";
+            node._node->Print(os); 
+          }
+          else if (node._intPoint && node._intPoint->_node)
+          {
+            os << "Node at intersection point: ";
+            node._intPoint->_node->Print(os); // intersection point
+          }
+          else
+            os << "mesh node is null\n";
+
+          return os;
+        }
+      };
+
+      std::vector< _nodeDef >              _nodes;
+      std::vector< int >                   _quantities;
+      _volumeDef*                          _next; // to store several _volumeDefs in a chain
+      StdMeshers::Cartesian3D::TGeomID     _solidID;
+      double                               _size;
+      const SMDS_MeshElement*              _volume; // new volume
+      std::vector<const SMDS_MeshElement*> _brotherVolume; // produced due to poly split
+      std::vector< SMESH_Block::TShapeID > _names; // name of side a polygon originates from
+
+      _volumeDef(): _next(0), _solidID(0), _size(0), _volume(0) {}
+      ~_volumeDef() { delete _next; }
+      _volumeDef( _volumeDef& other ):
+        _next(0), _solidID( other._solidID ), _size( other._size ), _volume( other._volume )
+      { _nodes.swap( other._nodes ); _quantities.swap( other._quantities ); other._volume = 0;
+        _names.swap( other._names ); }
+
+      size_t size() const { return 1 + ( _next ? _next->size() : 0 ); } // nb _volumeDef in a chain
+      _volumeDef* at(int index)
+      { return index == 0 ? this : ( _next ? _next->at(index-1) : _next ); }
+
+      void Set( _Node** nodes, int nb )
+      { _nodes.assign( nodes, nodes + nb ); }
+
+      void SetNext( _volumeDef* vd )
+      { if ( _next ) { _next->SetNext( vd ); } else { _next = vd; }}
+
+      bool IsEmpty() const { return (( _nodes.empty() ) &&
+                                     ( !_next || _next->IsEmpty() )); }
+      bool IsPolyhedron() const { return ( !_quantities.empty() ||
+                                           ( _next && !_next->_quantities.empty() )); }
+
+      std::vector<std::set<std::pair<int, int>>> getPolygonsEdges() const;
+      std::vector<std::set<std::pair<int, int>>> findOpenEdges() const;
+      int getStartNodeIndex(const int polygon) const;
+      std::map<int, std::vector<int>> findOverlappingPolygons() const;
+      bool divideOverlappingPolygons();
+      bool fixOpenEdgesPolygons();
+      bool capOpenEdgesPolygons(const std::vector<std::set<std::pair<int, int>>>& edgesByPolygon);
+      bool removeOpenEdgesPolygons(const std::vector<std::set<std::pair<int, int>>>& edgesByPolygon);
+
+      struct _linkDef: public std::pair<_ptr,_ptr> // to join polygons in removeExcessSideDivision()
+      {
+        _nodeDef _node1;//, _node2;
+        mutable /*const */_linkDef *_prev, *_next;
+        size_t _loopIndex;
+
+        _linkDef():_prev(0), _next(0) {}
+
+        void init( const _nodeDef& n1, const _nodeDef& n2, size_t iLoop )
+        {
+          _node1     = n1; //_node2 = n2;
+          _loopIndex = iLoop;
+          first      = n1.Ptr();
+          second     = n2.Ptr();
+          if ( first > second ) std::swap( first, second );
+        }
+        void setNext( _linkDef* next )
+        {
+          _next = next;
+          next->_prev = this;
+        }
+
+        friend std::ostream& operator<<(std::ostream& os, const _linkDef& link)
+        {
+          os << "Link def:\n";
+
+          os << link._node1;
+          if (link.first)
+            os << "first: " << link.first;
+
+          if (link.second)
+            os << "second: " << link.second;
+
+          os << "_loopIndex: " << link._loopIndex << '\n';
+
+          return os;
+        }
+      };
+    };
+
+    // topology of a hexahedron
+    static const std::size_t HEX_NODES_NUM = 8;
+    static const std::size_t HEX_LINKS_NUM = 12;
+    static const std::size_t HEX_QUADS_NUM = 6;
+    _Node _hexNodes [HEX_NODES_NUM];
+    _Link _hexLinks [HEX_LINKS_NUM];
+    _Face _hexQuads [HEX_QUADS_NUM];
+
+    // faces resulted from hexahedron intersection
+    std::vector< _Face > _polygons;
+
+    // intresections with EDGEs
+    std::vector< const StdMeshers::Cartesian3D::E_IntersectPoint* > _eIntPoints;
+
+    // additional nodes created at intersection points
+    std::vector< _Node > _intNodes;
+
+    // nodes inside the hexahedron (at VERTEXes) refer to _intNodes
+    std::vector< _Node* > _vIntNodes;
+
+    // computed volume elements
+    _volumeDef _volumeDefs;
+
+    StdMeshers::Cartesian3D::Grid*       _grid;
+    double      _sideLength[3];
+    int         _nbCornerNodes, _nbFaceIntNodes, _nbBndNodes;
+    int         _origNodeInd; // index of _hexNodes[0] node within the _grid
+    size_t      _i,_j,_k;
+    bool        _hasTooSmall;
+    int         _cellID;
+
+  public:
+    Hexahedron(StdMeshers::Cartesian3D::Grid* grid);
+    int MakeElements(SMESH_MesherHelper&                      helper,
+                     const TEdge2faceIDsMap& edge2faceIDsMap,
+                     const int numOfThreads = 1 );
+    void computeElements( const Solid* solid = 0, int solidIndex = -1 );
+
+  private:
+    Hexahedron(const Hexahedron& other, size_t i, size_t j, size_t k, int cellID );
+    void init( size_t i, size_t j, size_t k, const Solid* solid=0 );
+    void init( size_t i );
+    void clearNodesLinkedToNull(const Solid* solid, SMESH_MesherHelper& helper);
+    bool isSplittedLink(const Solid* solid, SMESH_MesherHelper& helper, const Hexahedron::_Link& linkIn) const;
+    void setIJK( size_t i );
+    /*Auxiliary methods to extract operations from monolitic compute method*/
+    void defineHexahedralFaces( const Solid* solid, const IsInternalFlag intFlag );
+    bool compute( const Solid* solid, const IsInternalFlag intFlag );
+    size_t getSolids( StdMeshers::Cartesian3D::TGeomID ids[] );
+    bool isCutByInternalFace( IsInternalFlag & maxFlag );
+    void addEdges(SMESH_MesherHelper&         helper,
+                  std::vector< Hexahedron* >& intersectedHex,
+                  const TEdge2faceIDsMap&     edge2faceIDsMap);
+    gp_Pnt findIntPoint( double u1, double proj1, double u2, double proj2,
+                         double proj, BRepAdaptor_Curve& curve,
+                         const gp_XYZ& axis, const gp_XYZ& origin );
+    int  getEntity( const StdMeshers::Cartesian3D::E_IntersectPoint* ip, int* facets, int& sub );
+    bool addIntersection( const StdMeshers::Cartesian3D::E_IntersectPoint* ip,
+                          std::vector< Hexahedron* >&  hexes,
+                          int ijk[], int dIJK[] );
+    bool isQuadOnFace( const size_t iQuad );
+    bool findChain( _Node* n1, _Node* n2, _Face& quad, std::vector<_Node*>& chainNodes );
+    bool closePolygon( _Face* polygon, std::vector<_Node*>& chainNodes ) const;
+    bool findChainOnEdge( const std::vector< _OrientedLink >& splits,
+                          const _OrientedLink&                prevSplit,
+                          const _OrientedLink&                avoidSplit,
+                          const std::set< StdMeshers::Cartesian3D::TGeomID > & concaveFaces,
+                          size_t &                            iS,
+                          _Face&                              quad,
+                          std::vector<_Node*>&                chn);
+    typedef std::pair< StdMeshers::Cartesian3D::TGeomID, int > TFaceOfLink; // (face, link)
+    static TFaceOfLink findStartLink(const std::vector< _OrientedLink* >& freeLinks,
+                                     std::set< StdMeshers::Cartesian3D::TGeomID >& usedFaceIDs);
+    size_t findCoplanarPolygon
+           (const _Face& thePolygon,
+            const size_t nbQuadPolygons,
+            std::vector< _OrientedLink* >& freeLinks,
+            int& nbFreeLinks,
+            const E_IntersectPoint& ipTmp,
+            std::set< StdMeshers::Cartesian3D::TGeomID >& usedFaceIDs,
+            std::map< StdMeshers::Cartesian3D::TGeomID, std::vector< const B_IntersectPoint* > >& tmpAddedFace,
+            const StdMeshers::Cartesian3D::TGeomID& curFace);
+    int  addVolumes( SMESH_MesherHelper& helper );
+    void addFaces( SMESH_MesherHelper& helper,
+                   const std::vector< const SMDS_MeshElement* > & boundaryVolumes );
+    void addSegments( SMESH_MesherHelper&     helper,
+                      const TEdge2faceIDsMap& edge2faceIDsMap );
+    void getVolumes( std::vector< const SMDS_MeshElement* > & volumes );
+    void getBoundaryElems( std::vector< const SMDS_MeshElement* > & boundaryVolumes );
+    void removeExcessSideDivision(const std::vector< Hexahedron* >& allHexa);
+    void removeExcessNodes(std::vector< Hexahedron* >& allHexa);
+    void preventVolumesOverlapping();
+    StdMeshers::Cartesian3D::TGeomID getAnyFace() const;
+    void cutByExtendedInternal( std::vector< Hexahedron* >& hexes,
+                                const TColStd_MapOfInteger& intEdgeIDs );
+    gp_Pnt mostDistantInternalPnt( int hexIndex, const gp_Pnt& p1, const gp_Pnt& p2 );
+    bool isOutPoint( _Link& link, int iP, SMESH_MesherHelper& helper, const Solid* solid ) const;
+    void sortVertexNodes(std::vector<_Node*>& nodes,
+                         _Node* curNode,
+                         StdMeshers::Cartesian3D::TGeomID face);
+    bool isInHole() const;
+    bool hasStrangeEdge() const;
+    bool checkPolyhedronSize( bool isCutByInternalFace, double & volSize ) const;
+    int checkPolyhedronValidity( _volumeDef* volDef, std::vector<std::vector<int>>& splitQuantities,
+                                 std::vector<std::vector<const SMDS_MeshNode*>>& splitNodes );
+    const SMDS_MeshElement* addPolyhedronToMesh( _volumeDef* volDef,
+                                                 SMESH_MesherHelper& helper,
+                                                 const std::vector<const SMDS_MeshNode*>& nodes,
+                                                 const std::vector<int>& quantities );
+    bool addHexa ();
+    bool addTetra();
+    bool addPenta();
+    bool addPyra ();
+    bool debugDumpLink( _Link* link );
+    _Node* findEqualNode( std::vector< _Node* >&       nodes,
+                          const StdMeshers::Cartesian3D::E_IntersectPoint* ip,
+                          const double            tol2 )
+    {
+      for ( size_t i = 0; i < nodes.size(); ++i )
+        if ( nodes[i]->EdgeIntPnt() == ip ||
+             nodes[i]->Point().SquareDistance( ip->_point ) <= tol2 )
+          return nodes[i];
+      return 0;
+    }
+    bool isCorner( const _Node* node ) const { return ( node >= &_hexNodes[0] &&
+                                                        node -  &_hexNodes[0] < 8 ); }
+    bool hasEdgesAround( const ConcaveFace* cf ) const;
+    bool isImplementEdges() const { return _grid->_edgeIntPool.nbElements(); }
+    bool isOutParam(const double uvw[3]) const;
+
+    typedef boost::container::flat_map< StdMeshers::Cartesian3D::TGeomID, size_t > TID2Nb;
+    static void insertAndIncrement( StdMeshers::Cartesian3D::TGeomID id, TID2Nb& id2nbMap )
+    {
+      TID2Nb::value_type s0( id, 0 );
+      TID2Nb::iterator id2nb = id2nbMap.insert( s0 ).first;
+      id2nb->second++;
+    }
+  }; // class Hexahedron
+} // end namespace Cartesian3D
+} // end namespace StdMeshers
+
+#endif

src/StdMeshers/StdMeshers_Import_1D2D.cxx

@@ -728,9 +728,11 @@ bool StdMeshers_Import_1D2D::Compute(SMESH_Mesh & theMesh, const TopoDS_Shape &
           {
             const SMDS_MeshNode* n = face->GetNode( i );
             if ( n == startNode || !checkedNodes.insert( n ).second ) continue;
-            helper.CheckNodeU( seamEdge, n, u=0, projTol, /*force=*/true );
-            nodesOnSeam.push_back( make_pair( u, n ));
-            ++nbNodesOnSeam;
+            if(helper.CheckNodeU( seamEdge, n, u=0, projTol, /*force=*/true ))
+            {
+              nodesOnSeam.push_back( make_pair( u, n ));
+              ++nbNodesOnSeam;
+            }
           }
         }
       }

src/StdMeshers/StdMeshers_Propagation.cxx

@@ -579,6 +579,14 @@ namespace {
     case HAS_PROPAG_HYP: {  // propag hyp on this submesh
       // --------------------------------------------------------
       switch ( event ) {
+      case SMESH_subMesh::ADD_FATHER_ALGO:
+      {
+        DBGMSG("HAS_PROPAG_HYP propagation to ADD_FATHER_ALGO " << subMesh->GetId());
+
+        // Rebuild propagation chain after an algo was added on father submesh
+        buildPropagationChain(subMesh);
+        break;
+      }
       case SMESH_subMesh::REMOVE_HYP:
       case SMESH_subMesh::REMOVE_FATHER_HYP: // remove propagation hyp
         if ( isPropagHyp && !getProagationHyp( subMesh ))

src/StdMeshers/StdMeshers_Regular_1D.cxx

@@ -492,6 +492,105 @@ static void compensateError(double a1, double an,
   }
 }
 
+
+//================================================================================
+/*!
+ * \brief adjust internal node parameters so that the last segment length == an,
+ *        and by distributing the error for the total length of curve segments
+ *        in relation to the target length computed from the current parameters
+ *  \param a1 - the first segment length
+ *  \param an - the last segment length
+ *  \param U1 - the first edge parameter
+ *  \param Un - the last edge parameter
+ *  \param length - the edge length
+ *  \param C3d - the edge curve
+ *  \param theParams - internal node parameters to adjust
+ */
+//================================================================================
+
+static void distributeError(double a1, double an,
+                            double U1, double Un,
+                            double            length,
+                            Adaptor3d_Curve&  C3d,
+                            list<double> &    theParams)
+{
+  // Compute the error of the total length based in the current curve parameters
+  double tol   = Min( Precision::Confusion(), 0.01 * Min(a1, an) );
+  double totalLength = 0.0;
+  double prevParam = U1;
+  list<double> segLengths;
+  list<double>::iterator itU = theParams.begin();
+  for ( ; itU != theParams.end(); ++itU )
+  {
+    // Compute the curve length between two adjacent parameters and sum them up
+    double curLength = GCPnts_AbscissaPoint::Length(C3d, prevParam, *itU, tol);
+    segLengths.push_back(curLength);
+    totalLength += curLength;
+    prevParam = *itU;
+  }
+  // Calculate the error between the total length of all segments based on given parameters
+  // and the target length of the edge itself
+  double error = totalLength - length;
+  // Compute the sum of all internal segments (= total computed length minus the length of
+  // the start and end segments)
+  double midLength = totalLength - (a1 + an);
+
+  // We only need to distribute the error, if the current parametrization is not correct,
+  // and if there are multiple internal segments
+  smIdType nPar = theParams.size();
+  if ( a1 + an <= length && nPar > 1 && fabs(error) > tol )
+  {
+    // Update the length of each internal segment (start and end length are given and not changed)
+    double newTotalLength = 0.0;
+    double newLength;
+    double relError = error / midLength;
+    list<double> newSegLengths;
+    list<double>::iterator itL = segLengths.begin();
+    for ( ; itL != segLengths.end(); ++itL )
+    {
+      // Do not update, but copy the first and the last segment lengths
+      newLength = *itL;
+      if (itL != segLengths.begin() && itL != --segLengths.end())
+      {
+        newLength -= newLength * relError;
+      }
+      newSegLengths.push_back(newLength);
+      newTotalLength += newLength;
+    }
+    bool reverse = ( U1 > Un );
+
+    // Update the parameters of the curve based on the new lengths
+    double curveLength, tol2, U;
+    double prevU = U1;
+    itU = theParams.begin();
+    itL = newSegLengths.begin();
+    for ( ; itU != theParams.end(); ++itU, ++itL )
+    {
+      curveLength = (reverse ? -(*itL) : *itL);
+      tol2        = Min( Precision::Confusion(), fabs(curveLength) / 100. );
+      GCPnts_AbscissaPoint Discret( tol2, C3d, curveLength, prevU );
+      if ( !Discret.IsDone() )
+      {
+        return;
+      }
+      U = Discret.Parameter();
+
+      double sign = reverse ? -1 : 1;
+      if ( sign*U1 < sign*U && sign*U < sign*Un )
+      {
+        *itU = U;
+      }
+      else
+      {
+        *itU = (sign*U >= sign*Un ? Un : U1);
+        break;
+      }
+      prevU = U;
+    }
+  }
+}
+
+
 //================================================================================
 /*!
  * \brief Class used to clean mesh on edges when 0D hyp modified.
@@ -1042,8 +1141,9 @@ bool StdMeshers_Regular_1D::computeInternalParameters(SMESH_Mesh &     theMesh,
       return error ( SMESH_Comment("Invalid segment lengths (")<<a1<<" and "<<an<<") "<<
                      "for an edge of length "<<theLength);
 
-    double q = ( an - a1 ) / ( 2 *theLength/( a1 + an ) - 1 );
-    int    n = int(fabs(q) > numeric_limits<double>::min() ? ( 1+( an-a1 )/q ) : ( 1+theLength/a1 ));
+    // Compute first the number of segments and then the arithmetic increment based on that number
+    int    n = static_cast<int>(2 * theLength / ( a1 + an ) + 0.5);
+    double q = (n > 1 ? ( an - a1 ) / (n - 1) : 0.0);
 
     double      U1 = theReverse ? l : f;
     double      Un = theReverse ? f : l;
@@ -1054,19 +1154,23 @@ bool StdMeshers_Regular_1D::computeInternalParameters(SMESH_Mesh &     theMesh,
       eltSize = -eltSize;
       q = -q;
     }
-    while ( n-- > 0 && eltSize * ( Un - U1 ) > 0 ) {
+    for (int i=0; i<n; i++) {
       // computes a point on a curve <theC3d> at the distance <eltSize>
       // from the point of parameter <param>.
       GCPnts_AbscissaPoint Discret( tol, theC3d, eltSize, param );
       if ( !Discret.IsDone() ) break;
       param = Discret.Parameter();
-      if ( param > f && param < l )
-        theParams.push_back( param );
-      else
-        break;
+      theParams.push_back( param );
       eltSize += q;
     }
-    compensateError( a1, an, U1, Un, theLength, theC3d, theParams );
+
+    distributeError( a1, an, U1, Un, theLength, theC3d, theParams );
+
+    // Do not include the parameter for the start or end of an edge in the list of parameters
+    // NOTE: it is required to correctly distribute the error
+    if (fabs(theParams.front() - U1) < tol) theParams.pop_front();
+    if (fabs(theParams.back() - Un) < tol)  theParams.pop_back();
+
     if ( theReverse ) theParams.reverse(); // NPAL18025
 
     return true;

src/StdMeshers/StdMeshers_ViscousLayers2D.cxx

@@ -1018,7 +1018,7 @@ bool _ViscousBuilder2D::makePolyLines()
     faceBndBox2D.Add( *_polyLineVec[ iPoLine]._segTree->getBox() );
   const double boxTol = 1e-3 * sqrt( faceBndBox2D.SquareExtent() );
 
-  if ( _maxThickness * maxLen2dTo3dRatio > sqrt( faceBndBox2D.SquareExtent() ) / 10 )
+  if ( _maxThickness * maxLen2dTo3dRatio > sqrt( faceBndBox2D.SquareExtent() ) / 2 )
   {
     vector< const _Segment* > foundSegs;
     double maxPossibleThick = 0;
@@ -1041,6 +1041,8 @@ bool _ViscousBuilder2D::makePolyLines()
           foundSegs.clear();
           L2._segTree->GetSegmentsNear( L1._lEdges[iLE]._ray, foundSegs );
           for ( size_t i = 0; i < foundSegs.size(); ++i )
+            // In this block periodically finding strange intersection with parameter close to 0
+            // to investigate
             if ( intersection.Compute( *foundSegs[i], L1._lEdges[iLE]._ray ))
             {
               double  distToL2 = intersection._param2 / L1._lEdges[iLE]._len2dTo3dRatio;
@@ -1865,26 +1867,30 @@ bool _ViscousBuilder2D::shrink()
         Geom2dAdaptor_Curve edgeCurve( pcurve, Min( uf, ul ), Max( uf, ul ));
         Geom2dAdaptor_Curve seg2Curve( seg2Line );
         Geom2dInt_GInter     curveInt( edgeCurve, seg2Curve, 1e-7, 1e-7 );
+
+        // In the older version length2D was set to this value only inside the !convex if block
+        // But it seems that length2D can be set here anyway, because if not set valid value of length2D here,
+        // it will be calculated later using length1D, and it can be not valid in cases if length1D is too large or too small.
+        length2D = L2->_lEdges[iFSeg2]._length2D;
+
+        /*  convex VERTEX
+         *                     L  seg2
+         *                    |  o---o---
+         *                    | /    |
+         *                    |/     |  L2
+         *                    x------x---      */
+        /* concave VERTEX
+         *                     o-----o---
+         *                      \    |
+         *                       \   |  L2
+         *                        x--x---
+         *                       /
+         *                    L /               */
         isConvex = ( curveInt.IsDone() && !curveInt.IsEmpty() );
-        if ( isConvex ) {
-          /*                   convex VERTEX */
+        if ( isConvex )
+        {
           length1D = Abs( u - curveInt.Point( 1 ).ParamOnFirst() );
-          double maxDist2d = 2 * L2->_lEdges[ iLSeg2 ]._length2D;
-          isConvex = ( length1D < maxDist2d * len1dTo2dRatio );
-          /*                                          |L  seg2
-           *                                          |  o---o---
-           *                                          | /    |
-           *                                          |/     |  L2
-           *                                          x------x---      */
-        }
-        if ( !isConvex ) { /* concave VERTEX */   /*  o-----o---
-                                                   *   \    |
-                                                   *    \   |  L2
-                                                   *     x--x---
-                                                   *    /
-                                                   * L /               */
-          length2D = L2->_lEdges[ iFSeg2 ]._length2D;
-          //if ( L2->_advancable ) continue;
+          length2D = Max(length2D, length1D / len1dTo2dRatio);
         }
       }
       else // L2 is advancable but in the face adjacent by L

src/Tools/smesh_plugins.py

@@ -61,6 +61,22 @@ except Exception as e:
   salome_pluginsmanager.logger.info('ERROR: MGCleaner plug-in is unavailable: {}'.format(e))
   pass
 
+try:
+  from mmgplugin.mmgPlug_plugin import Mmg
+  salome_pluginsmanager.AddFunction('ReMesh with MMG',
+                                    'Run MMG',
+                                    Mmg)
+except Exception as e:
+  salome_pluginsmanager.logger.info('ERROR: MMG plug-in is unavailable: {}'.format(e))
+  pass
+
+try:
+  from meshbooleanplugin.mesh_boolean_plugin import MeshBoolean
+  salome_pluginsmanager.AddFunction('Boolean Mesh Operations', 'Perform boolean operations on meshes', MeshBoolean)
+except Exception as e:
+  salome_pluginsmanager.logger.info('ERROR: MeshBoolean plug-in is unavailable: {}'.format(e))
+  pass
+
 try:
   from blocFissure.ihm.fissureCoude_plugin import fissureCoudeDlg
   salome_pluginsmanager.AddFunction('Meshed Pipe with a crack (blocFissure plugin)',
@@ -88,7 +104,7 @@ try:
                                       'Add a crack in a mesh with Zcracks plug-in',
                                       ZcracksLct)
 except Exception as e:
-  #print 'probleme zcracks'
+  #print 'problem zcracks'
   salome_pluginsmanager.logger.info('ERROR: Zcrack plug-in is unavailable: {}'.format(e))
   pass
 
@@ -99,6 +115,32 @@ try:
                                     'run topological volumic mesher',
                                     TopIIVolMeshLct)
 except Exception as e:
-  #print 'probleme zcracks'
   salome_pluginsmanager.logger.info('ERROR: TopIIVolMesh plug-in is unavailable: {}'.format(e))
   pass
+
+import sys
+import importlib
+import importlib.util
+
+smesh_pyplugin_dir = os.getenv("SMESH_PYPLUGIN_DIR")
+if smesh_pyplugin_dir:
+  plugin_dirs = smesh_pyplugin_dir.split(":")
+  # reverse plugin_dirs list to call smesh python plugin in good order
+  plugin_dirs.reverse()
+  for plg_dir in plugin_dirs:
+    if os.path.exists(plg_dir):
+      sys.path.insert(0, plg_dir)
+      #for filename in sorted(
+      #        filter(lambda x: os.path.isfile(os.path.join(plg_dir, x)),
+      #        os.listdir(plg_dir))):
+      for filename in os.listdir(plg_dir):
+        if filename.endswith(".py"):
+            f = os.path.join(plg_dir, filename)
+            try:
+              module_name = os.path.splitext(os.path.basename(f))[0]
+              _specs = importlib.util.find_spec(module_name)
+              _module = importlib.util.module_from_spec(_specs)
+              _specs.loader.exec_module(_module)
+              _module.init()
+            except:
+              print("Can not load python plugin from {}".format(plg_dir))

test/CTestTestfileInstall.cmake

@@ -41,3 +41,10 @@ FOREACH(tfile ${CPP_TESTS})
   ADD_TEST(${TEST_NAME} ${BASE_NAME} )
   SET_TESTS_PROPERTIES(${TEST_NAME} PROPERTIES LABELS "${COMPONENT_NAME};${COMPONENT_NAME}_tests")
 ENDFOREACH()
+
+FOREACH(tfile ${UNIT_TESTS})
+  GET_FILENAME_COMPONENT(BASE_NAME ${tfile} NAME_WE)
+  SET(TEST_NAME SMESH_${BASE_NAME})
+  ADD_TEST(${TEST_NAME} ${BASE_NAME} )
+  SET_TESTS_PROPERTIES(${TEST_NAME} PROPERTIES LABELS "${COMPONENT_NAME};${COMPONENT_NAME}_tests")
+ENDFOREACH()

test/SMESH_algo_switch_box.py

@@ -0,0 +1,81 @@
+# Tests that switching of algorithms back and forth does not lead to errors
+
+import salome
+salome.salome_init()
+
+from salome.geom import geomBuilder
+
+import  SMESH
+from salome.smesh import smeshBuilder
+
+# Create a box
+geompy = geomBuilder.New()
+
+O = geompy.MakeVertex(0, 0, 0)
+OX = geompy.MakeVectorDXDYDZ(1, 0, 0)
+OY = geompy.MakeVectorDXDYDZ(0, 1, 0)
+OZ = geompy.MakeVectorDXDYDZ(0, 0, 1)
+Box_1 = geompy.MakeBoxDXDYDZ(200, 200, 200)
+edge = geompy.CreateGroup(Box_1, geompy.ShapeType["EDGE"])
+geompy.UnionIDs(edge, [26])
+[edge] = geompy.GetExistingSubObjects(Box_1, False)
+geompy.addToStudy( O, 'O' )
+geompy.addToStudy( OX, 'OX' )
+geompy.addToStudy( OY, 'OY' )
+geompy.addToStudy( OZ, 'OZ' )
+geompy.addToStudy( Box_1, 'Box_1' )
+geompy.addToStudyInFather( Box_1, edge, 'edge' )
+
+# Create a mesh from the box and a sub-mesh from an edge
+smesh = smeshBuilder.New()
+Mesh_1 = smesh.Mesh(Box_1,'Mesh_1')
+Regular_1D = Mesh_1.Segment()
+Number_of_Segments_1 = Regular_1D.NumberOfSegments(15)
+Quadrangle_2D = Mesh_1.Quadrangle(algo=smeshBuilder.QUADRANGLE)
+Hexa_3D = Mesh_1.Hexahedron(algo=smeshBuilder.Hexa)
+edge_1 = Mesh_1.GroupOnGeom(edge,'edge',SMESH.EDGE)
+Regular_1D_1 = Mesh_1.Segment(geom=edge)
+Number_of_Segments_2 = Regular_1D_1.NumberOfSegments(2)
+Propagation_of_1D_Hyp = Regular_1D_1.Propagation()
+
+# Compute initial mesh
+Mesh_1.Compute()
+Mesh_1.CheckCompute()
+Sub_mesh_1 = Regular_1D_1.GetSubMesh()
+
+# Get the number of faces in the mesh
+num_faces_before = Mesh_1.NbFaces()
+print('Number of faces before switching: %d' % num_faces_before)
+
+# Switch to composite segment algorithm and compute the mesh
+status = Mesh_1.RemoveHypothesis(Regular_1D)
+CompositeSegment_1D = Mesh_1.Segment(algo=smeshBuilder.COMPOSITE)
+Mesh_1.AddHypothesis(CompositeSegment_1D)
+isDone = Mesh_1.Compute()
+Mesh_1.CheckCompute() # if propagation doesn't work it already fails here
+
+# Switch back to regular segment algorithm and compute the mesh
+status = Mesh_1.RemoveHypothesis(CompositeSegment_1D)
+Mesh_1.AddHypothesis(Regular_1D)
+Mesh_1.Compute()
+Mesh_1.CheckCompute()
+
+# Get the number of faces in the mesh
+num_faces_after = Mesh_1.NbFaces()
+print('Number of faces after switching: %d' % num_faces_after)
+assert num_faces_before == num_faces_after, 'Number of faces before and after switching should be the same'
+
+## Set names of Mesh objects
+smesh.SetName(CompositeSegment_1D.GetAlgorithm(), 'CompositeSegment_1D')
+smesh.SetName(Number_of_Segments_1, 'Number of Segments_1')
+smesh.SetName(Mesh_1.GetMesh(), 'Mesh_1')
+smesh.SetName(Number_of_Segments_2, 'Number of Segments_2')
+smesh.SetName(edge_1, 'edge')
+smesh.SetName(Hexa_3D.GetAlgorithm(), 'Hexa_3D')
+smesh.SetName(Sub_mesh_1, 'Sub-mesh_1')
+smesh.SetName(Regular_1D.GetAlgorithm(), 'Regular_1D')
+smesh.SetName(Propagation_of_1D_Hyp, 'Propagation of 1D Hyp. on Opposite Edges_1')
+smesh.SetName(Quadrangle_2D.GetAlgorithm(), 'Quadrangle_2D')
+
+if salome.sg.hasDesktop():
+  salome.sg.updateObjBrowser()

test/SMESH_algo_switch_face.py

@@ -0,0 +1,80 @@
+# Tests that switching of algorithms back and forth does not lead to errors
+
+import salome
+salome.salome_init()
+
+from salome.geom import geomBuilder
+
+import  SMESH
+from salome.smesh import smeshBuilder
+
+# Create a simple face
+geompy = geomBuilder.New()
+
+O = geompy.MakeVertex(0, 0, 0)
+OX = geompy.MakeVectorDXDYDZ(1, 0, 0)
+OY = geompy.MakeVectorDXDYDZ(0, 1, 0)
+OZ = geompy.MakeVectorDXDYDZ(0, 0, 1)
+Face_1 = geompy.MakeFaceHW(100, 100, 1)
+edge = geompy.CreateGroup(Face_1, geompy.ShapeType['EDGE'])
+geompy.UnionIDs(edge, [6])
+[edge] = geompy.GetExistingSubObjects(Face_1, False)
+geompy.addToStudy( O, 'O' )
+geompy.addToStudy( OX, 'OX' )
+geompy.addToStudy( OY, 'OY' )
+geompy.addToStudy( OZ, 'OZ' )
+geompy.addToStudy( Face_1, 'Face_1' )
+geompy.addToStudyInFather( Face_1, edge, 'edge' )
+
+# Create a mesh from the face and a sub-mesh from an edge
+smesh = smeshBuilder.New()
+
+Mesh_1 = smesh.Mesh(Face_1,'Mesh_1')
+Regular_1D = Mesh_1.Segment()
+Number_of_Segments_1 = Regular_1D.NumberOfSegments(3)
+Quadrangle_2D = Mesh_1.Quadrangle(algo=smeshBuilder.QUADRANGLE)
+edge_1 = Mesh_1.GroupOnGeom(edge,'edge',SMESH.EDGE)
+Regular_1D_1 = Mesh_1.Segment(geom=edge)
+Number_of_Segments_2 = Regular_1D_1.NumberOfSegments(2)
+Propagation_of_1D_Hyp = Regular_1D_1.Propagation()
+
+# Compute initial mesh
+Mesh_1.Compute()
+Mesh_1.CheckCompute()
+Sub_mesh_1 = Regular_1D_1.GetSubMesh()
+
+# Get the number of faces in the mesh
+num_faces_before = Mesh_1.NbFaces()
+print('Number of faces before switching: %d' % num_faces_before)
+
+# Switch to composite segment algorithm and compute the mesh
+status = Mesh_1.RemoveHypothesis(Regular_1D)
+CompositeSegment_1D = smesh.CreateHypothesis('CompositeSegment_1D')
+Mesh_1.AddHypothesis(CompositeSegment_1D)
+Mesh_1.Compute()
+Mesh_1.CheckCompute()
+
+# Switch back to regular segment algorithm and compute the mesh
+status = Mesh_1.RemoveHypothesis(CompositeSegment_1D)
+Mesh_1.AddHypothesis(Regular_1D)
+Mesh_1.Compute()
+Mesh_1.CheckCompute()
+
+# Get the number of faces in the mesh
+num_faces_after = Mesh_1.NbFaces()
+print('Number of faces after switching: %d' % num_faces_after)
+assert num_faces_before == num_faces_after, 'Number of faces before and after switching should be the same'
+
+## Set names of Mesh objects
+smesh.SetName(CompositeSegment_1D, 'CompositeSegment_1D')
+smesh.SetName(Number_of_Segments_1, 'Number of Segments_1')
+smesh.SetName(Mesh_1.GetMesh(), 'Mesh_1')
+smesh.SetName(Number_of_Segments_2, 'Number of Segments_2')
+smesh.SetName(edge_1, 'edge')
+smesh.SetName(Sub_mesh_1, 'Sub-mesh_1')
+smesh.SetName(Regular_1D.GetAlgorithm(), 'Regular_1D')
+smesh.SetName(Propagation_of_1D_Hyp, 'Propagation of 1D Hyp. on Opposite Edges_1')
+smesh.SetName(Quadrangle_2D.GetAlgorithm(), 'Quadrangle_2D')
+
+if salome.sg.hasDesktop():
+  salome.sg.updateObjBrowser()

test/data/HexahedronTest/NRTM1.brep

@@ -0,0 +1,441 @@
+DBRep_DrawableShape
+
+CASCADE Topology V1, (c) Matra-Datavision
+Locations 0
+Curve2ds 36
+1 0 0 1 0 
+1 0 0 1 0 
+1 10 0 0 -1 
+1 0 0 0 1 
+1 0 -10 1 0 
+1 0 0 1 0 
+1 0 0 0 -1 
+1 0 0 0 1 
+1 0 0 0 1 
+1 0 0 1 0 
+1 0 10 1 0 
+1 0 0 1 0 
+1 10 0 0 1 
+1 0 0 1 0 
+1 0 10 1 0 
+1 10 0 0 1 
+1 10 0 0 1 
+1 10 0 0 -1 
+1 0 0 0 1 
+1 0 10 1 0 
+1 0 10 1 0 
+1 0 -10 1 0 
+1 10 0 0 1 
+1 0 0 0 -1 
+2 5 5 1 0 -0 1 1
+1 0 6 1 0 
+2 5 5 1 0 -0 1 1
+1 0 6 1 0 
+1 6.2831853071795862 -0 0 1 
+1 0 -0 0 1 
+1 0 9 1 0 
+2 0 0 1 0 -0 1 1
+1 6.2831853071795862 -0 0 1 
+1 0 -0 0 1 
+1 0 0 1 0 
+2 0 0 1 0 -0 1 1
+Curves 17
+1 -5 -5 0 0 0 1 
+1 -5 -5 10 -0 1 0 
+1 -5 5 0 0 0 1 
+1 -5 -5 0 -0 1 0 
+1 -5 -5 0 1 0 -0 
+1 5 -5 0 0 0 1 
+1 -5 -5 10 1 0 -0 
+1 -5 5 10 1 0 -0 
+1 5 -5 10 -0 1 0 
+1 -5 5 0 1 0 -0 
+1 5 5 0 0 0 1 
+1 5 -5 0 -0 1 0 
+2 0 0 0 0 0 1 1 0 -0 -0 1 0 1
+1 1 -2.4492935982947064e-16 -6 0 0 1 
+2 0 0 3 0 0 1 1 0 -0 -0 1 0 1
+1 1 -2.4492935982947064e-16 -6 0 0 1 
+2 0 0 -6 0 0 1 1 0 -0 -0 1 0 1
+Polygon3D 0
+PolygonOnTriangulations 0
+Surfaces 11
+1 -5 -5 0 1 0 -0 0 0 1 0 -1 0 
+1 -5 -5 0 -0 1 0 0 0 1 1 0 -0 
+1 -5 -5 10 0 0 1 1 0 -0 -0 1 0 
+1 -5 5 0 -0 1 0 0 0 1 1 0 -0 
+1 -5 -5 0 0 0 1 1 0 -0 -0 1 0 
+1 5 -5 0 1 0 -0 0 0 1 0 -1 0 
+2 0 0 -6 0 0 1 1 0 -0 -0 1 0 1
+1 -5 -5 0 0 0 1 1 0 -0 -0 1 0 
+2 0 0 -6 0 0 1 1 0 -0 -0 1 0 1
+1 0 0 3 0 0 1 1 0 -0 -0 1 0 
+1 0 0 -6 0 0 1 1 0 -0 -0 1 0 
+Triangulations 0
+
+TShapes 58
+Ve
+1e-07
+-5 -5 10
+0 0
+
+0101101
+*
+Ve
+1e-07
+-5 -5 0
+0 0
+
+0101101
+*
+Ed
+ 1e-07 1 1 0
+1  1 0 0 10
+2  1 1 0 0 10
+2  2 2 0 0 10
+0
+
+0101000
+-58 0 +57 0 *
+Ve
+1e-07
+-5 5 10
+0 0
+
+0101101
+*
+Ed
+ 1e-07 1 1 0
+1  2 0 0 10
+2  3 1 0 0 10
+2  4 3 0 0 10
+0
+
+0101000
+-55 0 +58 0 *
+Ve
+1e-07
+-5 5 0
+0 0
+
+0101101
+*
+Ed
+ 1e-07 1 1 0
+1  3 0 0 10
+2  5 1 0 0 10
+2  6 4 0 0 10
+0
+
+0101000
+-55 0 +53 0 *
+Ed
+ 1e-07 1 1 0
+1  4 0 0 10
+2  7 1 0 0 10
+2  8 5 0 0 10
+0
+
+0101000
+-53 0 +57 0 *
+Wi
+
+0101100
+-56 0 -54 0 +52 0 +51 0 *
+Fa
+0  1e-07 1 0
+
+0101000
++50 0 *
+Ve
+1e-07
+5 -5 0
+0 0
+
+0101101
+*
+Ed
+ 1e-07 1 1 0
+1  5 0 0 10
+2  9 2 0 0 10
+2  10 5 0 0 10
+0
+
+0101000
+-48 0 +57 0 *
+Ve
+1e-07
+5 -5 10
+0 0
+
+0101101
+*
+Ed
+ 1e-07 1 1 0
+1  6 0 0 10
+2  11 2 0 0 10
+2  12 6 0 0 10
+0
+
+0101000
+-46 0 +48 0 *
+Ed
+ 1e-07 1 1 0
+1  7 0 0 10
+2  13 2 0 0 10
+2  14 3 0 0 10
+0
+
+0101000
+-46 0 +58 0 *
+Wi
+
+0101100
+-47 0 -45 0 +44 0 +56 0 *
+Fa
+0  1e-07 2 0
+
+0101000
++43 0 *
+Ve
+1e-07
+5 5 10
+0 0
+
+0101101
+*
+Ed
+ 1e-07 1 1 0
+1  8 0 0 10
+2  15 3 0 0 10
+2  16 4 0 0 10
+0
+
+0101000
+-41 0 +55 0 *
+Ed
+ 1e-07 1 1 0
+1  9 0 0 10
+2  17 3 0 0 10
+2  18 6 0 0 10
+0
+
+0101000
+-41 0 +46 0 *
+Wi
+
+0101100
+-54 0 -40 0 +39 0 +44 0 *
+Fa
+0  1e-07 3 0
+
+0101000
++38 0 *
+Ve
+1e-07
+5 5 0
+0 0
+
+0101101
+*
+Ed
+ 1e-07 1 1 0
+1  10 0 0 10
+2  19 4 0 0 10
+2  20 5 0 0 10
+0
+
+0101000
+-36 0 +53 0 *
+Ed
+ 1e-07 1 1 0
+1  11 0 0 10
+2  21 4 0 0 10
+2  22 6 0 0 10
+0
+
+0101000
+-41 0 +36 0 *
+Wi
+
+0101100
+-35 0 -34 0 +40 0 +52 0 *
+Fa
+0  1e-07 4 0
+
+0101000
++33 0 *
+Ed
+ 1e-07 1 1 0
+1  12 0 0 10
+2  23 5 0 0 10
+2  24 6 0 0 10
+0
+
+0101000
+-36 0 +48 0 *
+Wi
+
+0101100
+-51 0 -35 0 +47 0 +31 0 *
+Ve
+1.00000000244929e-07
+1 -2.44929359829471e-16 0
+0 0
+
+0101101
+*
+Ed
+ 1e-07 1 1 0
+1  13 0 0 6.28318530717959
+2  25 5 0 0 6.28318530717959
+2  26 7 0 0 6.28318530717959
+2  27 8 0 0 6.28318530717959
+2  28 9 0 0 6.28318530717959
+0
+
+0101000
++29 0 -29 0 *
+Wi
+
+0101100
+-28 0 *
+Fa
+0  1e-07 5 0
+
+0101000
++30 0 +27 0 *
+Wi
+
+0101100
+-45 0 -39 0 +34 0 +31 0 *
+Fa
+0  1e-07 6 0
+
+0101000
++25 0 *
+Ve
+1e-07
+1 -2.44929359829471e-16 3
+0 0
+
+0101101
+*
+Ed
+ 1e-07 1 1 0
+1  14 0 6 9
+3  29 30CN 7 0 6 9
+0
+
+0101000
++29 0 -23 0 *
+Ed
+ 1e-07 1 1 0
+1  15 0 0 6.28318530717959
+2  31 7 0 0 6.28318530717959
+2  32 10 0 0 6.28318530717959
+0
+
+0101000
++23 0 -23 0 *
+Wi
+
+0101100
++22 0 +28 0 -22 0 -21 0 *
+Fa
+0  1e-07 7 0
+
+0101000
++20 0 *
+Wi
+
+0101100
++21 0 *
+Fa
+0  1e-07 10 0
+
+0101000
++18 0 *
+Sh
+
+0101100
+-49 0 -42 0 +37 0 +32 0 -26 0 +24 0 -19 0 -17 0 *
+So
+
+0100000
++16 0 *
+Wi
+
+0101100
++28 0 *
+Fa
+0  1e-07 8 0
+
+0101000
++14 0 *
+Sh
+
+0101100
+-13 0 +19 0 +17 0 *
+So
+
+0100000
++12 0 *
+Ve
+1e-07
+1 -2.44929359829471e-16 -6
+0 0
+
+0101101
+*
+Ed
+ 1e-07 1 1 0
+1  16 0 0 6
+3  33 34CN 9 0 0 6
+0
+
+0101000
++10 0 -29 0 *
+Ed
+ 1e-07 1 1 0
+1  17 0 0 6.28318530717959
+2  35 9 0 0 6.28318530717959
+2  36 11 0 0 6.28318530717959
+0
+
+0101000
++10 0 -10 0 *
+Wi
+
+0101100
+-28 0 +9 0 +8 0 -9 0 *
+Fa
+0  1e-07 9 0
+
+0101000
++7 0 *
+Wi
+
+0101100
+-8 0 *
+Fa
+0  1e-07 11 0
+
+0101000
+-5 0 *
+Sh
+
+0101100
++6 0 +13 0 -4 0 *
+So
+
+0100000
++3 0 *
+Co
+
+1100000
++15 0 +11 0 +2 0 *
+
++1 0 
+0
+

test/data/HexahedronTest/NRTMJ4.brep

@@ -0,0 +1,598 @@
+DBRep_DrawableShape
+
+CASCADE Topology V1, (c) Matra-Datavision
+Locations 0
+Curve2ds 54
+1 0 0 1 0 
+1 0 0 1 0 
+1 100 0 0 -1 
+1 0 0 0 1 
+1 0 -50 1 0 
+1 0 0 1 0 
+1 0 0 0 -1 
+1 0 0 0 1 
+1 0 60 1 0 
+1 0 0 1 0 
+1 16.899999999999999 50 0 1 
+1 0 0 1 0 
+1 16.899999999999999 50 1 0 
+1 0 0 1 0 
+1 76.799999999999997 50 0 1 
+1 0 0 1 0 
+1 0 60 1 0 
+1 0 0 1 0 
+1 100 0 0 1 
+1 0 0 1 0 
+1 0 0 0 1 
+1 0 0 1 0 
+1 0 50 1 0 
+1 100 0 0 1 
+1 60 0 0 1 
+1 100 0 0 -1 
+1 0 60 1 0 
+1 0 -50 1 0 
+1 16.899999999999999 50 0 1 
+1 0 50 1 0 
+1 16.899999999999999 50 1 0 
+1 0 -50 1 0 
+1 76.799999999999997 50 0 1 
+1 0 50 1 0 
+1 0 60 1 0 
+1 0 -50 1 0 
+1 0 0 0 1 
+1 0 50 1 0 
+1 60 0 0 1 
+1 0 0 0 -1 
+2 25 25 -1 0 0 1 10
+1 0 0 1 0 
+1 16.899999999999999 0 0 -1 
+1 10 0 0 1 
+1 0 0 0 1 
+1 0 0 0 -1 
+1 59.899999999999999 0 0 -1 
+1 0 0 0 1 
+1 10 0 0 1 
+1 76.799999999999997 0 0 -1 
+1 0 50 1 0 
+2 0 0 1 0 -0 1 10
+1 6.2831853071795862 -0 0 1 
+1 0 -0 0 1 
+Curves 27
+1 0 0 0 0 0 1 
+1 0 0 100 -0 1 0 
+1 0 50 0 0 0 1 
+1 0 0 0 -0 1 0 
+1 60 0 0 0 0 1 
+1 50 0 16.899999999999999 1 0 -0 
+1 50 0 16.899999999999999 0 0 1 
+1 50 0 76.799999999999997 1 0 -0 
+1 60 0 0 0 0 1 
+1 0 0 100 1 0 -0 
+1 0 0 0 1 0 -0 
+1 0 50 100 1 0 -0 
+1 60 0 100 -0 1 0 
+1 60 50 0 0 0 1 
+1 50 50 16.899999999999999 1 0 -0 
+1 50 50 16.899999999999999 0 0 1 
+1 50 50 76.799999999999997 1 0 -0 
+1 60 50 0 0 0 1 
+1 0 50 0 1 0 -0 
+1 60 0 0 -0 1 0 
+2 25 25 0 0 0 -1 -1 0 -0 0 1 0 10
+1 60 0 16.899999999999999 -0 1 0 
+1 50 0 16.899999999999999 -0 1 0 
+1 50 0 76.799999999999997 -0 1 0 
+1 60 0 76.799999999999997 -0 1 0 
+2 25 25 -50 0 0 -1 -1 0 -0 0 1 0 10
+1 15.000000000000002 24.999999999999996 0 0 0 -1 
+Polygon3D 0
+PolygonOnTriangulations 0
+Surfaces 12
+1 0 0 0 1 0 -0 0 0 1 0 -1 0 
+1 0 0 0 -0 1 0 0 0 1 1 0 -0 
+1 0 0 100 0 0 1 1 0 -0 -0 1 0 
+1 0 50 0 -0 1 0 0 0 1 1 0 -0 
+1 0 0 0 0 0 1 1 0 -0 -0 1 0 
+1 60 0 0 1 0 -0 0 0 1 0 -1 0 
+1 50 0 16.899999999999999 0 0 1 1 0 -0 -0 1 0 
+1 50 0 16.899999999999999 1 0 -0 0 0 1 0 -1 0 
+1 50 0 76.799999999999997 0 0 1 1 0 -0 -0 1 0 
+1 60 0 0 1 0 -0 0 0 1 0 -1 0 
+2 25 25 0 0 0 -1 -1 0 -0 0 1 0 10
+1 25 25 -50 0 0 -1 -1 0 -0 0 1 0 
+Triangulations 0
+
+TShapes 72
+Ve
+1e-07
+0 0 100
+0 0
+
+0101101
+*
+Ve
+1e-07
+0 0 0
+0 0
+
+0101101
+*
+Ed
+ 1e-07 1 1 0
+1  1 0 0 100
+2  1 1 0 0 100
+2  2 2 0 0 100
+0
+
+0101000
+-72 0 +71 0 *
+Ve
+1e-07
+0 50 100
+0 0
+
+0101101
+*
+Ed
+ 1e-07 1 1 0
+1  2 0 0 50
+2  3 1 0 0 50
+2  4 3 0 0 50
+0
+
+0101000
+-69 0 +72 0 *
+Ve
+1e-07
+0 50 0
+0 0
+
+0101101
+*
+Ed
+ 1e-07 1 1 0
+1  3 0 0 100
+2  5 1 0 0 100
+2  6 4 0 0 100
+0
+
+0101000
+-69 0 +67 0 *
+Ed
+ 1e-07 1 1 0
+1  4 0 0 50
+2  7 1 0 0 50
+2  8 5 0 0 50
+0
+
+0101000
+-67 0 +71 0 *
+Wi
+
+0101100
+-70 0 -68 0 +66 0 +65 0 *
+Fa
+0  1e-07 1 0
+
+0101000
++64 0 *
+Ve
+1e-07
+60 0 0
+0 0
+
+0101101
+*
+Ve
+1e-07
+60 0 16.9
+0 0
+
+0101101
+*
+Ed
+ 1e-07 1 1 0
+1  5 0 0 16.9
+2  9 2 0 0 16.9
+2  10 6 0 0 16.9
+0
+
+0101000
++62 0 -61 0 *
+Ve
+1e-07
+50 0 16.9
+0 0
+
+0101101
+*
+Ed
+ 1e-07 1 1 0
+1  6 0 0 10
+2  11 2 0 0 10
+2  12 7 0 0 10
+0
+
+0101000
+-61 0 +59 0 *
+Ve
+1e-07
+50 0 76.8
+0 0
+
+0101101
+*
+Ed
+ 1e-07 1 1 0
+1  7 0 0 59.9
+2  13 2 0 0 59.9
+2  14 8 0 0 59.9
+0
+
+0101000
+-57 0 +59 0 *
+Ve
+1e-07
+60 0 76.8
+0 0
+
+0101101
+*
+Ed
+ 1e-07 1 1 0
+1  8 0 0 10
+2  15 2 0 0 10
+2  16 9 0 0 10
+0
+
+0101000
+-55 0 +57 0 *
+Ve
+1e-07
+60 0 100
+0 0
+
+0101101
+*
+Ed
+ 1e-07 1 1 0
+1  9 0 76.8 100
+2  17 2 0 76.8 100
+2  18 10 0 76.8 100
+0
+
+0101000
++55 0 -53 0 *
+Ed
+ 1e-07 1 1 0
+1  10 0 0 60
+2  19 2 0 0 60
+2  20 3 0 0 60
+0
+
+0101000
+-53 0 +72 0 *
+Ed
+ 1e-07 1 1 0
+1  11 0 0 60
+2  21 2 0 0 60
+2  22 5 0 0 60
+0
+
+0101000
+-62 0 +71 0 *
+Wi
+
+0101100
+-60 0 +58 0 -56 0 -54 0 -52 0 +51 0 +70 0 -50 0 *
+Fa
+0  1e-07 2 0
+
+0101000
++49 0 *
+Ve
+1e-07
+60 50 100
+0 0
+
+0101101
+*
+Ed
+ 1e-07 1 1 0
+1  12 0 0 60
+2  23 3 0 0 60
+2  24 4 0 0 60
+0
+
+0101000
+-47 0 +69 0 *
+Ed
+ 1e-07 1 1 0
+1  13 0 0 50
+2  25 3 0 0 50
+2  26 10 0 0 50
+0
+
+0101000
+-47 0 +53 0 *
+Wi
+
+0101100
+-68 0 -46 0 +45 0 +51 0 *
+Fa
+0  1e-07 3 0
+
+0101000
++44 0 *
+Ve
+1e-07
+60 50 0
+0 0
+
+0101101
+*
+Ve
+1e-07
+60 50 16.9
+0 0
+
+0101101
+*
+Ed
+ 1e-07 1 1 0
+1  14 0 0 16.9
+2  27 4 0 0 16.9
+2  28 6 0 0 16.9
+0
+
+0101000
++42 0 -41 0 *
+Ve
+1e-07
+50 50 16.9
+0 0
+
+0101101
+*
+Ed
+ 1e-07 1 1 0
+1  15 0 0 10
+2  29 4 0 0 10
+2  30 7 0 0 10
+0
+
+0101000
+-41 0 +39 0 *
+Ve
+1e-07
+50 50 76.8
+0 0
+
+0101101
+*
+Ed
+ 1e-07 1 1 0
+1  16 0 0 59.9
+2  31 4 0 0 59.9
+2  32 8 0 0 59.9
+0
+
+0101000
+-37 0 +39 0 *
+Ve
+1e-07
+60 50 76.8
+0 0
+
+0101101
+*
+Ed
+ 1e-07 1 1 0
+1  17 0 0 10
+2  33 4 0 0 10
+2  34 9 0 0 10
+0
+
+0101000
+-35 0 +37 0 *
+Ed
+ 1e-07 1 1 0
+1  18 0 76.8 100
+2  35 4 0 76.8 100
+2  36 10 0 76.8 100
+0
+
+0101000
++35 0 -47 0 *
+Ed
+ 1e-07 1 1 0
+1  19 0 0 60
+2  37 4 0 0 60
+2  38 5 0 0 60
+0
+
+0101000
+-42 0 +67 0 *
+Wi
+
+0101100
+-40 0 +38 0 -36 0 -34 0 -33 0 +46 0 +66 0 -32 0 *
+Fa
+0  1e-07 4 0
+
+0101000
++31 0 *
+Ed
+ 1e-07 1 1 0
+1  20 0 0 50
+2  39 5 0 0 50
+2  40 6 0 0 50
+0
+
+0101000
+-42 0 +62 0 *
+Wi
+
+0101100
+-65 0 -32 0 +50 0 +29 0 *
+Ve
+1e-07
+15 25 0
+0 0
+
+0101101
+*
+Ed
+ 1e-07 1 1 0
+1  21 0 0 6.28318530717959
+2  41 5 0 0 6.28318530717959
+2  42 11 0 0 6.28318530717959
+0
+
+0101000
++27 0 -27 0 *
+Wi
+
+0101100
++26 0 *
+Fa
+0  1e-07 5 0
+
+0101000
++28 0 +25 0 *
+Ed
+ 1e-07 1 1 0
+1  22 0 0 50
+2  43 6 0 0 50
+2  44 7 0 0 50
+0
+
+0101000
+-41 0 +61 0 *
+Wi
+
+0101100
+-60 0 -23 0 +40 0 +29 0 *
+Fa
+0  1e-07 6 0
+
+0101000
++22 0 *
+Ed
+ 1e-07 1 1 0
+1  23 0 0 50
+2  45 7 0 0 50
+2  46 8 0 0 50
+0
+
+0101000
+-39 0 +59 0 *
+Wi
+
+0101100
+-20 0 -38 0 +23 0 +58 0 *
+Fa
+0  1e-07 7 0
+
+0101000
++19 0 *
+Ed
+ 1e-07 1 1 0
+1  24 0 0 50
+2  47 8 0 0 50
+2  48 9 0 0 50
+0
+
+0101000
+-37 0 +57 0 *
+Wi
+
+0101100
+-56 0 -17 0 +36 0 +20 0 *
+Fa
+0  1e-07 8 0
+
+0101000
++16 0 *
+Ed
+ 1e-07 1 1 0
+1  25 0 0 50
+2  49 9 0 0 50
+2  50 10 0 0 50
+0
+
+0101000
+-35 0 +55 0 *
+Wi
+
+0101100
+-17 0 -34 0 +14 0 +54 0 *
+Fa
+0  1e-07 9 0
+
+0101000
++13 0 *
+Wi
+
+0101100
+-52 0 -45 0 +33 0 +14 0 *
+Fa
+0  1e-07 10 0
+
+0101000
++11 0 *
+Ve
+1e-07
+15 25 -50
+0 0
+
+0101101
+*
+Ed
+ 1e-07 1 1 0
+1  26 0 0 6.28318530717959
+2  51 11 0 0 6.28318530717959
+2  52 12 0 0 6.28318530717959
+0
+
+0101000
++9 0 -9 0 *
+Ed
+ 1e-07 1 1 0
+1  27 0 0 50
+3  53 54CN 11 0 0 50
+0
+
+0101000
+-9 0 +27 0 *
+Wi
+
+0101100
+-8 0 +7 0 +26 0 -7 0 *
+Fa
+0  1e-07 11 0
+
+0101000
++6 0 *
+Wi
+
+0101100
++8 0 *
+Fa
+0  1e-07 12 0
+
+0101000
++4 0 *
+Sh
+
+0101100
+-63 0 -48 0 +43 0 +30 0 -24 0 +21 0 +18 0 +15 0 -12 0 +10 0 
++5 0 +3 0 *
+So
+
+1100000
++2 0 *
+
++1 0 
+0
+

test/tests.set

@@ -109,6 +109,8 @@ SET(GOOD_TESTS
   ex31_dimGroup.py
   PAL_MESH_043_2D.py
   SMESH_AdvancedEditor.py
+  SMESH_algo_switch_box.py
+  SMESH_algo_switch_face.py
   SMESH_blocks.py
   SMESH_box.py
   SMESH_BuildCompound.py
@@ -127,6 +129,11 @@ SET(CPP_TESTS
   SMESH_RegularGridTest
 )
 
+SET(UNIT_TESTS # Any unit test add in src names space should be added here 
+  HexahedronTest
+  HexahedronCanonicalShapesTest
+)
+
 # The following tests can be executed without driver, just by python.
 # ----------------------------------------------------------------------------