Salome HOME
bos #20256: [CEA 18523] Porting SMESH to int 64 bits
[plugins/blsurfplugin.git] / src / BLSURFPlugin / BLSURFPlugin_BLSURF.cxx
index d593f6408f087a984e9d9059a707ebce459a5233..6e30bcaa83681fb42d766085af953dcbc419bea3 100644 (file)
@@ -1,4 +1,4 @@
-// Copyright (C) 2007-2016  CEA/DEN, EDF R&D
+// Copyright (C) 2007-2021  CEA/DEN, EDF R&D
 //
 // This library is free software; you can redistribute it and/or
 // modify it under the terms of the GNU Lesser General Public
@@ -58,6 +58,7 @@ extern "C"{
 #include <cstdlib>
 
 // OPENCASCADE includes
+#include <BRepBndLib.hxx>
 #include <BRepBuilderAPI_MakeFace.hxx>
 #include <BRepBuilderAPI_MakePolygon.hxx>
 #include <BRepBuilderAPI_MakeWire.hxx>
@@ -133,13 +134,13 @@ namespace
   static PyMethodDef PyStdOut_methods[] = {
     {"write",  (PyCFunction)PyStdOut_write,  METH_VARARGS,
       PyDoc_STR("write(string) -> None")},
-    {NULL,    NULL}   /* sentinel */
+    {0, 0, 0, 0}   /* sentinel */
   };
 
   static PyMemberDef PyStdOut_memberlist[] = {
     {(char*)"softspace", T_INT,  offsetof(PyStdOut, softspace), 0,
      (char*)"flag indicating that a space needs to be printed; used by print"},
-    {NULL} /* Sentinel */
+    {0, 0, 0, 0, 0} /* Sentinel */
   };
 
   static PyTypeObject PyStdOut_Type = {
@@ -187,6 +188,14 @@ namespace
     0,                            /*tp_new*/
     0,                            /*tp_free*/
     0,                            /*tp_is_gc*/
+    0,                            /*tp_bases*/
+    0,                            /*tp_mro*/
+    0,                            /*tp_cache*/
+    0,                            /*tp_subclasses*/
+    0,                            /*tp_weaklist*/
+    0,                            /*tp_del*/
+    0,                            /*tp_version_tag*/
+    0,                            /*tp_finalize*/
   };
 
   PyObject * newPyStdOut( std::string& out )
@@ -495,7 +504,7 @@ TopoDS_Shape BLSURFPlugin_BLSURF::entryToShape(std::string entry)
 {
   GEOM::GEOM_Object_var aGeomObj;
   TopoDS_Shape S = TopoDS_Shape();
-  SALOMEDS::SObject_var aSObj = SMESH_Gen_i::getStudyServant()->FindObjectID( entry.c_str() );
+  SALOMEDS::SObject_var aSObj = SMESH_Gen_i::GetSMESHGen()->getStudyServant()->FindObjectID( entry.c_str() );
   if (!aSObj->_is_nil()) {
     CORBA::Object_var obj = aSObj->GetObject();
     aGeomObj = GEOM::GEOM_Object::_narrow(obj);
@@ -612,7 +621,7 @@ void BLSURFPlugin_BLSURF::createEnforcedVertexOnFace(TopoDS_Shape faceShape, BLS
 /////////////////////////////////////////////////////////
 void createAttractorOnFace(TopoDS_Shape GeomShape, std::string AttractorFunction, double defaultSize)
 {
-  double xa, ya, za; // Coordinates of attractor point
+  double xa=0., ya=0., za=0.; // Coordinates of attractor point
   double a, b;       // Attractor parameter
   double d = 0.;
   bool createNode=false; // To create a node on attractor projection
@@ -879,6 +888,12 @@ void BLSURFPlugin_BLSURF::SetParameters(const BLSURFPlugin_Hypothesis* hyp,
   bool   _quadraticMesh         = BLSURFPlugin_Hypothesis::GetDefaultQuadraticMesh();
   int    _verb                  = BLSURFPlugin_Hypothesis::GetDefaultVerbosity();
   //int    _topology              = BLSURFPlugin_Hypothesis::GetDefaultTopology();
+  bool   _useSurfaceProximity      = BLSURFPlugin_Hypothesis::GetDefaultUseSurfaceProximity     ();
+  int    _nbSurfaceProximityLayers = BLSURFPlugin_Hypothesis::GetDefaultNbSurfaceProximityLayers();
+  double _surfaceProximityRatio    = BLSURFPlugin_Hypothesis::GetDefaultSurfaceProximityRatio   ();
+  bool   _useVolumeProximity       = BLSURFPlugin_Hypothesis::GetDefaultUseVolumeProximity      ();
+  int    _nbVolumeProximityLayers  = BLSURFPlugin_Hypothesis::GetDefaultNbVolumeProximityLayers ();
+  double _volumeProximityRatio     = BLSURFPlugin_Hypothesis::GetDefaultVolumeProximityRatio    ();
 
   // PreCAD
   //int _precadMergeEdges         = BLSURFPlugin_Hypothesis::GetDefaultPreCADMergeEdges();
@@ -941,6 +956,13 @@ void BLSURFPlugin_BLSURF::SetParameters(const BLSURFPlugin_Hypothesis* hyp,
     //_precadRemoveDuplicateCADFaces = hyp->GetPreCADRemoveDuplicateCADFaces();
     //_precadProcess3DTopology = hyp->GetPreCADProcess3DTopology();
     //_precadDiscardInput      = hyp->GetPreCADDiscardInput();
+    _useSurfaceProximity      = hyp->GetUseSurfaceProximity     ();
+    _nbSurfaceProximityLayers = hyp->GetNbSurfaceProximityLayers();
+    _surfaceProximityRatio    = hyp->GetSurfaceProximityRatio   ();
+    _useVolumeProximity       = hyp->GetUseVolumeProximity      ();
+    _nbVolumeProximityLayers  = hyp->GetNbVolumeProximityLayers ();
+    _volumeProximityRatio     = hyp->GetVolumeProximityRatio    ();
+
 
     const BLSURFPlugin_Hypothesis::TOptionValues& opts = hyp->GetOptionValues();
     BLSURFPlugin_Hypothesis::TOptionValues::const_iterator opIt;
@@ -969,6 +991,14 @@ void BLSURFPlugin_BLSURF::SetParameters(const BLSURFPlugin_Hypothesis* hyp,
       if ( !opIt->second.empty() ) {
         set_param(css, opIt->first.c_str(), opIt->second.c_str());
       }
+
+    if ( hyp->GetHyperPatches().size() < hyp->GetHyperPatchEntries().size() )
+    {
+      std::map< std::string, TopoDS_Shape > entryToShape;
+      FillEntryToShape( hyp, entryToShape );
+      const_cast<BLSURFPlugin_Hypothesis*>( hyp )->SetHyperPatchIDsByEntry( theGeomShape,
+                                                                            entryToShape );
+    }
   }
 
   if ( BLSURFPlugin_Hypothesis::HasPreCADOptions( hyp ))
@@ -988,11 +1018,12 @@ void BLSURFPlugin_BLSURF::SetParameters(const BLSURFPlugin_Hypothesis* hyp,
      case BLSURFPlugin_Hypothesis::PhysicalGlobalSize:
        set_param(css, "physical_size_mode", "global");
        set_param(css, "global_physical_size", _phySizeRel ? val_to_string_rel(_phySize).c_str() : val_to_string(_phySize).c_str());
+       //useGradation = true;
        break;
      case BLSURFPlugin_Hypothesis::PhysicalLocalSize:
        set_param(css, "physical_size_mode", "local");
        set_param(css, "global_physical_size", _phySizeRel ? val_to_string_rel(_phySize).c_str() : val_to_string(_phySize).c_str());
-       useGradation = true;
+       //useGradation = true;
        break;
      default:
        set_param(css, "physical_size_mode", "none");
@@ -1004,13 +1035,13 @@ void BLSURFPlugin_BLSURF::SetParameters(const BLSURFPlugin_Hypothesis* hyp,
        set_param(css, "geometric_size_mode", "global");
        set_param(css, "geometric_approximation", val_to_string(_angleMesh).c_str());
        set_param(css, "chordal_error", val_to_string(_chordalError).c_str());
-       useGradation = true;
+       //useGradation = true;
        break;
      case BLSURFPlugin_Hypothesis::GeometricalLocalSize:
        set_param(css, "geometric_size_mode", "local");
        set_param(css, "geometric_approximation", val_to_string(_angleMesh).c_str());
        set_param(css, "chordal_error", val_to_string(_chordalError).c_str());
-       useGradation = true;
+       //useGradation = true;
        break;
      default:
        set_param(css, "geometric_size_mode", "none");
@@ -1041,6 +1072,7 @@ void BLSURFPlugin_BLSURF::SetParameters(const BLSURFPlugin_Hypothesis* hyp,
      // - if maxsize is not explicitly specified, we pass default value computed automatically, in this case "relative" flag is ignored
      set_param(css, "max_size", _maxSizeRel ? val_to_string_rel(_maxSize).c_str() : val_to_string(_maxSize).c_str());
    }
+   useGradation = true; // bos #18758
    // anisotropic and quadrangle mesh requires disabling gradation
    if ( _anisotropic && _elementType != BLSURFPlugin_Hypothesis::Triangles )
      useGradation = false; // limitation of V1.3
@@ -1086,6 +1118,19 @@ void BLSURFPlugin_BLSURF::SetParameters(const BLSURFPlugin_Hypothesis* hyp,
    set_param(css, "element_order",                     _quadraticMesh ? "quadratic" : "linear");
    set_param(css, "verbose",                           val_to_string(_verb).c_str());
 
+   set_param(css, "use_surface_proximity",             _useSurfaceProximity ? "yes" : "no" );
+   if ( _useSurfaceProximity )
+   {
+     set_param(css, "surface_proximity_layers",        SMESH_Comment( _nbSurfaceProximityLayers ));
+     set_param(css, "surface_proximity_ratio",         SMESH_Comment( _surfaceProximityRatio ));
+   }
+   set_param(css, "use_volume_proximity",             _useVolumeProximity ? "yes" : "no" );
+   if ( _useVolumeProximity )
+   {
+     set_param(css, "volume_proximity_layers",        SMESH_Comment( _nbVolumeProximityLayers ));
+     set_param(css, "volume_proximity_ratio",         SMESH_Comment( _volumeProximityRatio ));
+   }
+
    _smp_phy_size = _phySizeRel ? _phySize*diagonal : _phySize;
    if ( _verb > 0 )
      std::cout << "_smp_phy_size = " << _smp_phy_size << std::endl;
@@ -1520,7 +1565,7 @@ namespace
       //double tol = (( u2node.rbegin()->first - u2node.begin()->first ) / 20.) / u2node.size();
       Standard_Real f,l;
       BRep_Tool::Range( TopoDS::Edge( shape ), f,l );
-      double tol = (( l - f ) / 20.) / u2node.size();
+      double tol = (( l - f ) / 10.) / double( u2node.size() ); // 10. - adjusted for #17262
 
       std::multimap< double, const SMDS_MeshNode* >::iterator un2, un1;
       for ( un2 = u2node.begin(), un1 = un2++; un2 != u2node.end(); un1 = un2++ )
@@ -1773,7 +1818,7 @@ namespace
           if ( !n2nIt->second ) {
             n->GetXYZ( xyz );
             gp_XY uv = tmpHelper.GetNodeUV( _proxyFace, n );
-            n2nIt->second = helper.AddNode( xyz[0], xyz[1], xyz[2], uv.X(), uv.Y() );
+            n2nIt->second = helper.AddNode( xyz[0], xyz[1], xyz[2], /*id=*/0, uv.X(), uv.Y() );
           }
           nodes[ nbN ] = n2nIt->second;
         }
@@ -1895,6 +1940,7 @@ bool BLSURFPlugin_BLSURF::compute(SMESH_Mesh&         aMesh,
   typedef set< SMESHDS_SubMesh*, ShapeTypeCompare > TSubMeshSet;
   TSubMeshSet edgeSubmeshes;
   TSubMeshSet& mergeSubmeshes = edgeSubmeshes;
+  double existingPhySize = 0;
 
   TopTools_IndexedMapOfShape pmap, emap, fmap;
 
@@ -2193,6 +2239,8 @@ bool BLSURFPlugin_BLSURF::compute(SMESH_Mesh&         aMesh,
           // tmin and tmax can change in case of viscous layer on an adjacent edge
           tmin = nodeDataVec.front().param;
           tmax = nodeDataVec.back().param;
+
+          existingPhySize += nodeData->Length() / double( nodeDataVec.size() - 1 );
         }
         else
         {
@@ -2240,7 +2288,7 @@ bool BLSURFPlugin_BLSURF::compute(SMESH_Mesh&         aMesh,
       if ( nodeData )
       {
         const std::vector<UVPtStruct>& nodeDataVec = nodeData->GetUVPtStruct();
-        const int                      nbNodes     = nodeDataVec.size();
+        const int                      nbNodes     = (int) nodeDataVec.size();
 
         dcad_edge_discretization_t *dedge;
         dcad_get_edge_discretization(dcad, edg, &dedge);
@@ -2342,44 +2390,6 @@ bool BLSURFPlugin_BLSURF::compute(SMESH_Mesh&         aMesh,
     } // for edge
   } //for face
 
-  // Clear mesh from already meshed edges if possible else
-  // remember that merge is needed
-  TSubMeshSet::iterator smIt = edgeSubmeshes.begin();
-  for ( ; smIt != edgeSubmeshes.end(); ++smIt ) // loop on already meshed EDGEs
-  {
-    SMESHDS_SubMesh* smDS = *smIt;
-    if ( !smDS ) continue;
-    SMDS_NodeIteratorPtr nIt = smDS->GetNodes();
-    if ( nIt->more() )
-    {
-      const SMDS_MeshNode* n = nIt->next();
-      if ( n->NbInverseElements( SMDSAbs_Face ) > 0 )
-      {
-        needMerge = true; // to correctly sew with viscous mesh
-        // add existing medium nodes to helper
-        if ( aMesh.NbEdges( ORDER_QUADRATIC ) > 0 )
-        {
-          SMDS_ElemIteratorPtr edgeIt = smDS->GetElements();
-          while ( edgeIt->more() )
-            helper.AddTLinks( static_cast<const SMDS_MeshEdge*>(edgeIt->next()));
-        }
-        continue;
-      }
-    }
-    if ( allowSubMeshClearing )
-    {
-      SMDS_ElemIteratorPtr eIt = smDS->GetElements();
-      while ( eIt->more() ) meshDS->RemoveFreeElement( eIt->next(), 0 );
-      SMDS_NodeIteratorPtr nIt = smDS->GetNodes();
-      while ( nIt->more() ) meshDS->RemoveFreeNode( nIt->next(), 0 );
-      smDS->Clear();
-    }
-    else
-    {
-      needMerge = true;
-    }
-  }
-
   ///////////////////////
   // PERIODICITY       //
   ///////////////////////
@@ -2405,8 +2415,9 @@ bool BLSURFPlugin_BLSURF::compute(SMESH_Mesh&         aMesh,
       {
         // If no source points, call periodicity without transformation function
         meshgems_cad_periodicity_transformation_t periodicity_transformation = NULL;
-        status = cad_add_face_multiple_periodicity_with_transformation_function(c, theFace1_ids_c, theFace1_ids.size(),
-                                                                                theFace2_ids_c, theFace2_ids.size(), periodicity_transformation, NULL);
+        status = cad_add_face_multiple_periodicity_with_transformation_function
+          (c, theFace1_ids_c, (meshgems_integer) theFace1_ids.size(), theFace2_ids_c,
+           (meshgems_integer) theFace2_ids.size(), periodicity_transformation, NULL);
         if(status != STATUS_OK)
           cout << "cad_add_face_multiple_periodicity_with_transformation_function failed with error code " << status << "\n";
       }
@@ -2415,11 +2426,12 @@ bool BLSURFPlugin_BLSURF::compute(SMESH_Mesh&         aMesh,
         // get the transformation vertices
         double* theSourceVerticesCoords_c = &_preCadFacesIDsPeriodicityVector[i].theSourceVerticesCoords[0];
         double* theTargetVerticesCoords_c = &_preCadFacesIDsPeriodicityVector[i].theTargetVerticesCoords[0];
-        int nbSourceVertices = _preCadFacesIDsPeriodicityVector[i].theSourceVerticesCoords.size()/3;
-        int nbTargetVertices = _preCadFacesIDsPeriodicityVector[i].theTargetVerticesCoords.size()/3;
+        int nbSourceVertices = (int)_preCadFacesIDsPeriodicityVector[i].theSourceVerticesCoords.size()/3;
+        int nbTargetVertices = (int)_preCadFacesIDsPeriodicityVector[i].theTargetVerticesCoords.size()/3;
 
-        status = cad_add_face_multiple_periodicity_with_transformation_function_by_points(c, theFace1_ids_c, theFace1_ids.size(),
-                                                                                          theFace2_ids_c, theFace2_ids.size(), theSourceVerticesCoords_c, nbSourceVertices, theTargetVerticesCoords_c, nbTargetVertices);
+        status = cad_add_face_multiple_periodicity_with_transformation_function_by_points
+          (c, theFace1_ids_c, (meshgems_integer) theFace1_ids.size(), theFace2_ids_c,
+           (meshgems_integer) theFace2_ids.size(), theSourceVerticesCoords_c, nbSourceVertices, theTargetVerticesCoords_c, nbTargetVertices);
         if(status != STATUS_OK)
           cout << "cad_add_face_multiple_periodicity_with_transformation_function_by_points failed with error code " << status << "\n";
       }
@@ -2450,8 +2462,9 @@ bool BLSURFPlugin_BLSURF::compute(SMESH_Mesh&         aMesh,
       {
         // If no source points, call periodicity without transformation function
         meshgems_cad_periodicity_transformation_t periodicity_transformation = NULL;
-        status = cad_add_edge_multiple_periodicity_with_transformation_function(c, theEdge1_ids_c, theEdge1_ids.size(),
-                                                                                theEdge2_ids_c, theEdge2_ids.size(), periodicity_transformation, NULL);
+        status = cad_add_edge_multiple_periodicity_with_transformation_function
+          (c, theEdge1_ids_c, (meshgems_integer) theEdge1_ids.size(), theEdge2_ids_c,
+           (meshgems_integer) theEdge2_ids.size(), periodicity_transformation, NULL);
         if(status != STATUS_OK)
           cout << "cad_add_edge_multiple_periodicity_with_transformation_function failed with error code " << status << "\n";
       }
@@ -2460,17 +2473,34 @@ bool BLSURFPlugin_BLSURF::compute(SMESH_Mesh&         aMesh,
         // get the transformation vertices
         double* theSourceVerticesCoords_c = &_preCadEdgesIDsPeriodicityVector[i].theSourceVerticesCoords[0];
         double* theTargetVerticesCoords_c = &_preCadEdgesIDsPeriodicityVector[i].theTargetVerticesCoords[0];
-        int nbSourceVertices = _preCadEdgesIDsPeriodicityVector[i].theSourceVerticesCoords.size()/3;
-        int nbTargetVertices = _preCadEdgesIDsPeriodicityVector[i].theTargetVerticesCoords.size()/3;
+        int nbSourceVertices = (int) _preCadEdgesIDsPeriodicityVector[i].theSourceVerticesCoords.size()/3;
+        int nbTargetVertices = (int) _preCadEdgesIDsPeriodicityVector[i].theTargetVerticesCoords.size()/3;
 
-        status = cad_add_edge_multiple_periodicity_with_transformation_function_by_points(c, theEdge1_ids_c, theEdge1_ids.size(),
-                                                                                          theEdge2_ids_c, theEdge2_ids.size(), theSourceVerticesCoords_c, nbSourceVertices, theTargetVerticesCoords_c, nbTargetVertices);
+        status = cad_add_edge_multiple_periodicity_with_transformation_function_by_points
+          (c, theEdge1_ids_c, (meshgems_integer) theEdge1_ids.size(), theEdge2_ids_c,
+           (meshgems_integer) theEdge2_ids.size(), theSourceVerticesCoords_c, nbSourceVertices, theTargetVerticesCoords_c, nbTargetVertices);
         if(status != STATUS_OK)
           cout << "cad_add_edge_multiple_periodicity_with_transformation_function_by_points failed with error code " << status << "\n";
       }
     }
   }
 
+  if ( !_hypothesis && !edgeSubmeshes.empty() && existingPhySize != 0 )
+  {
+    // prevent failure due to the default PhySize incompatible with size of existing 1D mesh
+    // and with face size
+    // double minFaceSize = existingPhySize / edgeSubmeshes.size();
+    // for ( int iF = 1; iF <= fmap.Extent(); ++iF )
+    // {
+    //   Bnd_Box box;
+    //   BRepBndLib::Add( fmap( iF ), box );
+    //   gp_XYZ delta = box.CornerMax().XYZ() - box.CornerMin().XYZ();
+    //   std::sort( delta.ChangeData(), delta.ChangeData() + 3 );
+    //   minFaceSize = Min( minFaceSize, delta.Coord(2) );
+    // }
+    // set_param(css, "global_physical_size", val_to_string( minFaceSize * 0.5 ).c_str());
+    // set_param(css, "max_size",             val_to_string( minFaceSize * 5 ).c_str());
+  }
   
   // TODO: be able to use a mesh in input.
   // See imsh usage in Products/templates/mg-cadsurf_template_common.cpp
@@ -2513,12 +2543,51 @@ bool BLSURFPlugin_BLSURF::compute(SMESH_Mesh&         aMesh,
 
   mesh_t *msh = NULL;
   cadsurf_get_mesh(css, &msh);
-  if(!msh){
+  if ( !msh || STATUS_IS_ERROR( status ))
+  {
     /* release the mesh object */
     cadsurf_regain_mesh(css, msh);
     return error(_comment);
   }
 
+  // Clear mesh from already meshed edges if possible else
+  // remember that merge is needed
+  TSubMeshSet::iterator smIt = edgeSubmeshes.begin();
+  for ( ; smIt != edgeSubmeshes.end(); ++smIt ) // loop on already meshed EDGEs
+  {
+    SMESHDS_SubMesh* smDS = *smIt;
+    if ( !smDS ) continue;
+    SMDS_NodeIteratorPtr nIt = smDS->GetNodes();
+    if ( nIt->more() )
+    {
+      const SMDS_MeshNode* n = nIt->next();
+      if ( n->NbInverseElements( SMDSAbs_Face ) > 0 )
+      {
+        needMerge = true; // to correctly sew with viscous mesh
+        // add existing medium nodes to helper
+        if ( aMesh.NbEdges( ORDER_QUADRATIC ) > 0 )
+        {
+          SMDS_ElemIteratorPtr edgeIt = smDS->GetElements();
+          while ( edgeIt->more() )
+            helper.AddTLinks( static_cast<const SMDS_MeshEdge*>(edgeIt->next()));
+        }
+        continue;
+      }
+    }
+    if ( allowSubMeshClearing )
+    {
+      SMDS_ElemIteratorPtr eIt = smDS->GetElements();
+      while ( eIt->more() ) meshDS->RemoveFreeElement( eIt->next(), 0 );
+      SMDS_NodeIteratorPtr nIt = smDS->GetNodes();
+      while ( nIt->more() ) meshDS->RemoveFreeNode( nIt->next(), 0 );
+      smDS->Clear();
+    }
+    else
+    {
+      needMerge = true;
+    }
+  }
+
   std::string GMFFileName = BLSURFPlugin_Hypothesis::GetDefaultGMFFile();
   if (_hypothesis)
     GMFFileName = _hypothesis->GetGMFFile();
@@ -2599,8 +2668,7 @@ bool BLSURFPlugin_BLSURF::compute(SMESH_Mesh&         aMesh,
           }
           if (!groupDone)
           {
-            int groupId;
-            SMESH_Group* aGroup = aMesh.AddGroup(SMDSAbs_Node, currentEnfVertex->grpName.c_str(), groupId);
+            SMESH_Group* aGroup = aMesh.AddGroup( SMDSAbs_Node, currentEnfVertex->grpName.c_str() );
             aGroup->SetName( currentEnfVertex->grpName.c_str() );
             SMESHDS_Group* aGroupDS = static_cast<SMESHDS_Group*>( aGroup->GetGroupDS() );
             aGroupDS->SMDSGroup().Add(nodes[iv]);
@@ -2701,6 +2769,8 @@ bool BLSURFPlugin_BLSURF::compute(SMESH_Mesh&         aMesh,
                             nodes[evtri[0]], nodes[evtri[1]], nodes[evtri[2]]);
     }
     else {
+      if ( helper.GetIsQuadratic() )
+        helper.SetSubShape( tag );
       tri = helper.AddFace(nodes[vtx[0]], nodes[vtx[1]], nodes[vtx[2]]);
     }
     meshDS->SetMeshElementOnShape(tri, tag);
@@ -2730,7 +2800,6 @@ bool BLSURFPlugin_BLSURF::compute(SMESH_Mesh&         aMesh,
     };
     if (type == MESHGEMS_MESH_ELEMENT_TYPE_QUAD9) {
       // QUADRATIC QUADRANGLE
-      std::cout << "This is a quadratic quadrangle" << std::endl;
       if (tags[evquad[0]]) {
         meshDS->SetNodeOnFace(nodes[evquad[0]], tag);
         tags[evquad[0]] = false;
@@ -2895,6 +2964,10 @@ bool BLSURFPlugin_BLSURF::Compute(SMESH_Mesh & aMesh, SMESH_MesherHelper* aHelpe
   context_t *ctx = context_new();
   if (!ctx) return error("Pb in context_new()");
 
+  if ( aMesh.NbNodes() > std::numeric_limits< meshgems_integer >::max() ||
+       aMesh.NbFaces() > std::numeric_limits< meshgems_integer >::max() )
+    return error("Too large input mesh");
+
   BLSURF_Cleaner cleaner( ctx );
 
   message_cb_user_data mcud;
@@ -2913,7 +2986,7 @@ bool BLSURFPlugin_BLSURF::Compute(SMESH_Mesh & aMesh, SMESH_MesherHelper* aHelpe
   // Fill an input mesh
 
   mesh_t * msh = meshgems_mesh_new_in_memory( ctx );
-  if ( !msh ) return error("Pb. in meshgems_mesh_new_in_memory()"); 
+  if ( !msh ) return error("Pb. in meshgems_mesh_new_in_memory()");
 
   // mark nodes used by 2D elements
   SMESHDS_Mesh* meshDS = aMesh.GetMeshDS();
@@ -2923,7 +2996,7 @@ bool BLSURFPlugin_BLSURF::Compute(SMESH_Mesh & aMesh, SMESH_MesherHelper* aHelpe
     const SMDS_MeshNode* n = nodeIt->next();
     n->setIsMarked( n->NbInverseElements( SMDSAbs_Face ));
   }
-  meshgems_mesh_set_vertex_count( msh, meshDS->NbNodes() );
+  meshgems_mesh_set_vertex_count( msh, (meshgems_integer) meshDS->NbNodes() );
 
   // set node coordinates
   if ( meshDS->NbNodes() != meshDS->MaxNodeID() )
@@ -2940,8 +3013,8 @@ bool BLSURFPlugin_BLSURF::Compute(SMESH_Mesh & aMesh, SMESH_MesherHelper* aHelpe
   }
 
   // set nodes of faces
-  meshgems_mesh_set_triangle_count  ( msh, meshDS->GetMeshInfo().NbTriangles() );
-  meshgems_mesh_set_quadrangle_count( msh, meshDS->GetMeshInfo().NbQuadrangles() );
+  meshgems_mesh_set_triangle_count  ( msh, (meshgems_integer) meshDS->GetMeshInfo().NbTriangles() );
+  meshgems_mesh_set_quadrangle_count( msh, (meshgems_integer) meshDS->GetMeshInfo().NbQuadrangles() );
   meshgems_integer nodeIDs[4];
   meshgems_integer iT = 1, iQ = 1;
   SMDS_FaceIteratorPtr faceIt = meshDS->facesIterator();
@@ -2952,7 +3025,7 @@ bool BLSURFPlugin_BLSURF::Compute(SMESH_Mesh & aMesh, SMESH_MesherHelper* aHelpe
     if ( nbNodes > 4 || face->IsPoly() ) continue;
 
     for ( i = 0; i < nbNodes; ++i )
-      nodeIDs[i] = face->GetNode( i )->GetID();
+      nodeIDs[i] = (meshgems_integer) face->GetNode( i )->GetID();
     if ( nbNodes == 3 )
       meshgems_mesh_set_triangle_vertices  ( msh, iT++, nodeIDs );
     else
@@ -3012,14 +3085,14 @@ bool BLSURFPlugin_BLSURF::Compute(SMESH_Mesh & aMesh, SMESH_MesherHelper* aHelpe
   {
     meshgems_mesh_get_vertex_coordinates( omsh, i, xyz );
     SMDS_MeshNode* n = meshDS->AddNode( xyz[0], xyz[1], xyz[2] );
-    nodeID = n->GetID();
+    nodeID = (meshgems_integer) n->GetID();
     meshgems_mesh_set_vertex_tag( omsh, i, &nodeID ); // save mapping of IDs in MG and SALOME meshes
   }
 
   // add triangles
   meshgems_integer nbtri = 0;
   meshgems_mesh_get_triangle_count( omsh, &nbtri );
-  const SMDS_MeshNode* nodes[3];
+  const SMDS_MeshNode* nodes[4];
   for ( i = 1; i <= nbtri; ++i )
   {
     meshgems_mesh_get_triangle_vertices( omsh, i, nodeIDs );
@@ -3031,6 +3104,33 @@ bool BLSURFPlugin_BLSURF::Compute(SMESH_Mesh & aMesh, SMESH_MesherHelper* aHelpe
     meshDS->AddFace( nodes[0], nodes[1], nodes[2] );
   }
 
+  // add quadrangles
+  meshgems_integer nbquad = 0;
+  meshgems_mesh_get_quadrangle_count( omsh, &nbquad );
+  for ( i = 1; i <= nbquad; ++i )
+  {
+    meshgems_mesh_get_quadrangle_vertices( omsh, i, nodeIDs );
+    for ( int j = 0; j < 4; ++j )
+    {
+      meshgems_mesh_get_vertex_tag( omsh, nodeIDs[j], &nodeID );
+      nodes[j] = meshDS->FindNode( nodeID );
+    }
+    meshDS->AddFace( nodes[0], nodes[1], nodes[2], nodes[3] );
+  }
+
+  if ( _hypothesis )
+  {
+    std::string GMFFileName = _hypothesis->GetGMFFile();
+    if ( !GMFFileName.empty() )
+    {
+      bool asciiFound  = (GMFFileName.find(".mesh", GMFFileName.size()-5) != std::string::npos);
+      bool binaryFound = (GMFFileName.find(".meshb",GMFFileName.size()-6) != std::string::npos);
+      if ( !asciiFound && !binaryFound )
+        GMFFileName.append(".mesh");
+      mesh_write_mesh(msh, GMFFileName.c_str());
+    }
+  }
+
   cadsurf_regain_mesh(css, omsh);
 
   // as we don't assign the new triangles to a shape (the pseudo-shape),
@@ -3325,7 +3425,7 @@ status_t message_cb(message_t *msg, void *user_data)
        err.find("periodicity") != string::npos )
   {
     // remove ^A from the tail
-    int len = strlen( desc );
+    size_t len = strlen( desc );
     while (len > 0 && desc[len-1] != '\n')
       len--;
     mcud->_error->append( desc, len );
@@ -3431,8 +3531,8 @@ bool BLSURFPlugin_BLSURF::Evaluate(SMESH_Mesh&         aMesh,
       nb1d = (int)( fullAng/_angleMesh + 1 );
     }
     fullNbSeg += nb1d;
-    std::vector<int> aVec(SMDSEntity_Last);
-    for(int i=SMDSEntity_Node; i<SMDSEntity_Last; i++) aVec[i]=0;
+    std::vector<smIdType> aVec(SMDSEntity_Last);
+    for(smIdType i=SMDSEntity_Node; i<SMDSEntity_Last; i++) aVec[i]=0;
     if( IsQuadratic > 0 ) {
       aVec[SMDSEntity_Node] = 2*nb1d - 1;
       aVec[SMDSEntity_Quad_Edge] = nb1d;
@@ -3479,7 +3579,7 @@ bool BLSURFPlugin_BLSURF::Evaluate(SMESH_Mesh&         aMesh,
         nbTria = nbQuad = nbTria / 3 + 1;
       }
     }
-    std::vector<int> aVec(SMDSEntity_Last,0);
+    std::vector<smIdType> aVec(SMDSEntity_Last,0);
     if( IsQuadratic ) {
       int nb1d_in = (nbTria*3 - nb1d) / 2;
       aVec[SMDSEntity_Node] = nbNodes + nb1d_in;
@@ -3503,8 +3603,8 @@ bool BLSURFPlugin_BLSURF::Evaluate(SMESH_Mesh&         aMesh,
   double tetrVol = 0.1179*ELen*ELen*ELen;
   int nbVols  = int(aVolume/tetrVol);
   int nb1d_in = int(( nbVols*6 - fullNbSeg ) / 6 );
-  std::vector<int> aVec(SMDSEntity_Last);
-  for(int i=SMDSEntity_Node; i<SMDSEntity_Last; i++) aVec[i]=0;
+  std::vector<smIdType> aVec(SMDSEntity_Last);
+  for(smIdType i=SMDSEntity_Node; i<SMDSEntity_Last; i++) aVec[i]=0;
   if( IsQuadratic ) {
     aVec[SMDSEntity_Node] = nb1d_in/3 + 1 + nb1d_in;
     aVec[SMDSEntity_Quad_Tetra] = nbVols;
@@ -3518,3 +3618,23 @@ bool BLSURFPlugin_BLSURF::Evaluate(SMESH_Mesh&         aMesh,
 
   return true;
 }
+
+//================================================================================
+/*!
+ * \brief Find TopoDS_Shape for each hyper-patch study entry in a hypothesis
+ */
+//================================================================================
+
+void BLSURFPlugin_BLSURF::FillEntryToShape( const BLSURFPlugin_Hypothesis*          hyp,
+                                            std::map< std::string, TopoDS_Shape > & entryToShape )
+{
+  SMESH_Gen_i* smeshGen = SMESH_Gen_i::GetSMESHGen();
+  for ( const ::BLSURFPlugin_Hypothesis::THyperPatchEntries& entries : hyp->GetHyperPatchEntries() )
+    for ( const std::string& entry : entries )
+    {
+      GEOM::GEOM_Object_var go = smeshGen->GetGeomObjectByEntry( entry );
+      TopoDS_Shape       shape = smeshGen->GeomObjectToShape( go );
+      if ( !shape.IsNull() )
+        entryToShape.insert({ entry, shape });
+    }
+}