Salome HOME
[bos #38045] [EDF] (2023-T3) Stand alone version for Netgen meshers.
authorcconopoima <cesar.conopoima@gmail.com>
Tue, 24 Oct 2023 15:58:28 +0000 (16:58 +0100)
committercesarconopoima <cesar.conopoima@opencascade.com>
Thu, 18 Jan 2024 10:30:12 +0000 (10:30 +0000)
Separating mesh face and proxy mesh for quand adaptors and VL in 3D.

Intermedial commit

First version of 1D SA.

2D SA version and definition of simpleHypo.

Define the NETGENPlugin_NETGEN_1D2D3D_SA class to implement 1D,1D2D,1D2D3D SA version of the mesh. Using NETGENPlugin_NETGEN_2D_SA class to define the 2D SA version.

Cleaning code.

Toward 2D SA version.

Function to feed 1D elements from imported mesh into netgen.

2D SA version for netgen.

Finish support for parallel computing of 2D mesh with netgen. Add the NETGEN_2D_Remote. UsIn NETGEN_2D_SA class use bounding box of faces to allow fast discarting edge nodes far from the faces been meshed.

23 files changed:
idl/NETGENPlugin_Algorithm.idl
src/NETGENPlugin/CMakeLists.txt
src/NETGENPlugin/NETGENPluginBuilder.py
src/NETGENPlugin/NETGENPlugin_DriverParam.cxx
src/NETGENPlugin/NETGENPlugin_DriverParam.hxx
src/NETGENPlugin/NETGENPlugin_Mesher.cxx
src/NETGENPlugin/NETGENPlugin_Mesher.hxx
src/NETGENPlugin/NETGENPlugin_NETGEN_1D2D3D_SA.cxx [new file with mode: 0644]
src/NETGENPlugin/NETGENPlugin_NETGEN_1D2D3D_SA.hxx [new file with mode: 0644]
src/NETGENPlugin/NETGENPlugin_NETGEN_2D_ONLY.cxx
src/NETGENPlugin/NETGENPlugin_NETGEN_2D_ONLY.hxx
src/NETGENPlugin/NETGENPlugin_NETGEN_2D_Remote.cxx [new file with mode: 0644]
src/NETGENPlugin/NETGENPlugin_NETGEN_2D_Remote.hxx [new file with mode: 0644]
src/NETGENPlugin/NETGENPlugin_NETGEN_2D_Remote_i.cxx [new file with mode: 0644]
src/NETGENPlugin/NETGENPlugin_NETGEN_2D_Remote_i.hxx [new file with mode: 0644]
src/NETGENPlugin/NETGENPlugin_NETGEN_2D_SA.cxx [new file with mode: 0644]
src/NETGENPlugin/NETGENPlugin_NETGEN_2D_SA.hxx [new file with mode: 0644]
src/NETGENPlugin/NETGENPlugin_NETGEN_3D.cxx
src/NETGENPlugin/NETGENPlugin_NETGEN_3D_Remote.cxx
src/NETGENPlugin/NETGENPlugin_NETGEN_3D_SA.cxx
src/NETGENPlugin/NETGENPlugin_Runner_main.cxx
src/NETGENPlugin/NETGENPlugin_SimpleHypothesis_2D.cxx
src/NETGENPlugin/NETGENPlugin_i.cxx

index 3b9499c64c38b76e8b32c5baefd0aba0e83385c9..10510e984ef6136740a7fc875e10341ea8c19fce 100644 (file)
@@ -74,6 +74,13 @@ module NETGENPlugin
   {
   };
 
+  /*!
+   * NETGENPlugin_NETGEN_2D: interface of "Remote Triangle (Netgen)" algorithm
+   */
+  interface NETGENPlugin_NETGEN_2D_Remote : NETGENPlugin::NETGENPlugin_NETGEN_2D_ONLY
+  {
+  };
+
   /*!
    * NETGENPlugin_Remesher_2D: interface of "NETGEN Remesher" algorithm,
    * generating 2D elements basing on an existing 2D mesh
index ee09f36152a15248bd2d49adfaa27273b8113fd1..b22d617c39b905c17e66b1aeff6304f731e983f6 100644 (file)
@@ -27,15 +27,16 @@ INCLUDE_DIRECTORIES(
   ${GUI_INCLUDE_DIRS}
   ${GEOM_INCLUDE_DIRS}
   ${OpenCASCADE_INCLUDE_DIR}
-  ${NETGEN_INCLUDE_DIRS}
   ${SMESH_INCLUDE_DIRS}
   ${MEDCOUPLING_INCLUDE_DIRS}
   ${Boost_INCLUDE_DIRS}
   ${OMNIORB_INCLUDE_DIR}
   ${PROJECT_BINARY_DIR}/idl
-  ${NETGEN_ZLIB_INCLUDE_DIRS}
 )
 
+#Avoid compilation warnings from netgen headers
+INCLUDE_DIRECTORIES( SYSTEM  ${NETGEN_INCLUDE_DIRS} ${NETGEN_ZLIB_INCLUDE_DIRS} )
+
 # additional preprocessor / compiler flags
 ADD_DEFINITIONS(
   ${QT_DEFINITIONS}
@@ -100,8 +101,12 @@ SET(NETGENEngine_HEADERS
   NETGENPlugin_Defs.hxx
   NETGENPlugin_DriverParam.hxx
   NETGENPlugin_NETGEN_3D_SA.hxx
+  NETGENPlugin_NETGEN_2D_SA.hxx
+  NETGENPlugin_NETGEN_1D2D3D_SA.hxx
   NETGENPlugin_NETGEN_3D_Remote.hxx
   NETGENPlugin_NETGEN_3D_Remote_i.hxx
+  NETGENPlugin_NETGEN_2D_Remote.hxx
+  NETGENPlugin_NETGEN_2D_Remote_i.hxx
 )
 
 # --- sources ---
@@ -131,8 +136,12 @@ SET(NETGENEngine_SOURCES
   NETGENPlugin_i.cxx
   NETGENPlugin_DriverParam.cxx
   NETGENPlugin_NETGEN_3D_SA.cxx
+  NETGENPlugin_NETGEN_2D_SA.cxx
+  NETGENPlugin_NETGEN_1D2D3D_SA.cxx
   NETGENPlugin_NETGEN_3D_Remote.cxx
   NETGENPlugin_NETGEN_3D_Remote_i.cxx
+  NETGENPlugin_NETGEN_2D_Remote.cxx
+  NETGENPlugin_NETGEN_2D_Remote_i.cxx
 )
 
 SET(NetgenRunner_SOURCES
index 7ee865527a89e364e7f7e93d3e8b205101326f12..9ff05ad58762e0640282e5f26fd0f72595b8e00a 100644 (file)
@@ -49,6 +49,7 @@ NETGEN_1D2D3D = "NETGEN_2D3D"
 NETGEN_1D2D   = "NETGEN_2D"
 ## Algorithm type: Netgen triangle 2D algorithm, see NETGEN_2D_Only_Algorithm
 NETGEN_2D     = "NETGEN_2D_ONLY"
+NETGEN_2D_Remote = "NETGEN_2D_Remote"
 ## Algorithm type: Synonim of NETGEN_1D2D3D, see NETGEN_1D2D3D_Algorithm
 NETGEN_FULL   = NETGEN_1D2D3D
 ## Algorithm type: Synonim of NETGEN_3D, see NETGEN_3D_Algorithm
@@ -510,6 +511,33 @@ class NETGEN_3D_Remote_Algorithm(NETGEN_3D_Algorithm):
 
     pass # end of NETGEN_3D_Remote_Algorithm class
 
+## Tetrahedron 2D algorithm
+#
+#  It can be created by calling smeshBuilder.Mesh.Triangle() or smeshBuilder.Mesh.Triangle( smeshBuilder.NETGEN, geom=0 )
+#
+class NETGEN_2D_Remote_Algorithm(NETGEN_2D_Only_Algorithm):
+
+    ## type of algorithm used with helper function in smeshBuilder.Mesh class
+    #  @internal
+    algoType   = NETGEN_2D_Remote
+    ## flag pointing either this algorithm should be used by default in dynamic method
+    #  of smeshBuilder.Mesh class
+    #  @internal
+    isDefault  = False
+    ## doc string of the method
+    #  @internal
+    docHelper  = "Remotely Creates triangles in face of solids"
+
+    ## Private constructor.
+    #  @param mesh parent mesh object algorithm is assigned to
+    #  @param geom geometry (shape/sub-shape) algorithm is assigned to;
+    #              if it is @c 0 (default), the algorithm is assigned to the main shape
+    def __init__(self, mesh, geom=0):
+        self.algoType = NETGEN_2D_Remote
+        NETGEN_2D_Only_Algorithm.__init__(self, mesh, geom)
+        pass
+
+    pass # end of NETGEN_2D_Remote_Algorithm class
 
 
 ## Triangle (helper) 1D-2D algorithm
index fedea933a748d4200e28720ab9a249c0a8383216..aa8502ac8f02bc8816a6630d3b98c28b0155212d 100644 (file)
  * @param aParams The object to display
  */
 void printNetgenParams(netgen_params& aParams){
-  std::cout << "has_netgen_param: " << aParams.has_netgen_param << std::endl;
-  std::cout << "maxh: " << aParams.maxh << std::endl;
-  std::cout << "minh: " << aParams.minh << std::endl;
-  std::cout << "segmentsperedge: " << aParams.segmentsperedge << std::endl;
-  std::cout << "grading: " << aParams.grading << std::endl;
-  std::cout << "curvaturesafety: " << aParams.curvaturesafety << std::endl;
-  std::cout << "secondorder: " << aParams.secondorder << std::endl;
-  std::cout << "quad: " << aParams.quad << std::endl;
-  std::cout << "optimize: " << aParams.optimize << std::endl;
-  std::cout << "fineness: " << aParams.fineness << std::endl;
-  std::cout << "uselocalh: " << aParams.uselocalh << std::endl;
-  std::cout << "merge_solids: " << aParams.merge_solids << std::endl;
-  std::cout << "chordalError: " << aParams.chordalError << std::endl;
-  std::cout << "optsteps2d: " << aParams.optsteps2d << std::endl;
-  std::cout << "optsteps3d: " << aParams.optsteps3d << std::endl;
-  std::cout << "elsizeweight: " << aParams.elsizeweight << std::endl;
-  std::cout << "opterrpow: " << aParams.opterrpow << std::endl;
-  std::cout << "delaunay: " << aParams.delaunay << std::endl;
-  std::cout << "checkoverlap: " << aParams.checkoverlap << std::endl;
-  std::cout << "checkchartboundary: " << aParams.checkchartboundary << std::endl;
-  std::cout << "closeedgefac: " << aParams.closeedgefac << std::endl;
-  std::cout << "nbThreadMesher: " << aParams.nbThreads << std::endl;
-  std::cout << "has_local_size: " << aParams.has_local_size << std::endl;
-  std::cout << "meshsizefilename: " << aParams.meshsizefilename << std::endl;
-  std::cout << "has_maxelementvolume_hyp: " << aParams.has_maxelementvolume_hyp << std::endl;
-  std::cout << "maxElementVolume: " << aParams.maxElementVolume << std::endl;
-  std::cout << "has_LengthFromEdges_hyp: " << aParams.has_LengthFromEdges_hyp << std::endl;
+  if ( aParams.myType == Hypo )
+  {
+    std::cout << "has_netgen_param: " << aParams.has_netgen_param << std::endl;
+    std::cout << "maxh: " << aParams.maxh << std::endl;
+    std::cout << "minh: " << aParams.minh << std::endl;
+    std::cout << "segmentsperedge: " << aParams.segmentsperedge << std::endl;
+    std::cout << "grading: " << aParams.grading << std::endl;
+    std::cout << "curvaturesafety: " << aParams.curvaturesafety << std::endl;
+    std::cout << "secondorder: " << aParams.secondorder << std::endl;
+    std::cout << "quad: " << aParams.quad << std::endl;
+    std::cout << "optimize: " << aParams.optimize << std::endl;
+    std::cout << "fineness: " << aParams.fineness << std::endl;
+    std::cout << "uselocalh: " << aParams.uselocalh << std::endl;
+    std::cout << "merge_solids: " << aParams.merge_solids << std::endl;
+    std::cout << "chordalError: " << aParams.chordalError << std::endl;
+    std::cout << "optsteps2d: " << aParams.optsteps2d << std::endl;
+    std::cout << "optsteps3d: " << aParams.optsteps3d << std::endl;
+    std::cout << "elsizeweight: " << aParams.elsizeweight << std::endl;
+    std::cout << "opterrpow: " << aParams.opterrpow << std::endl;
+    std::cout << "delaunay: " << aParams.delaunay << std::endl;
+    std::cout << "checkoverlap: " << aParams.checkoverlap << std::endl;
+    std::cout << "checkchartboundary: " << aParams.checkchartboundary << std::endl;
+    std::cout << "closeedgefac: " << aParams.closeedgefac << std::endl;
+    std::cout << "nbThreadMesher: " << aParams.nbThreads << std::endl;
+    std::cout << "has_local_size: " << aParams.has_local_size << std::endl;
+    std::cout << "meshsizefilename: " << aParams.meshsizefilename << std::endl;
+    std::cout << "has_maxelementvolume_hyp: " << aParams.has_maxelementvolume_hyp << std::endl;
+    std::cout << "maxElementVolume: " << aParams.maxElementVolume << std::endl;
+    std::cout << "has_LengthFromEdges_hyp: " << aParams.has_LengthFromEdges_hyp << std::endl;
+  }
 }
 
-/**
- * @brief Import a param_file into a netgen_params structure
- *
- * @param param_file Name of the file
- * @param aParams Structure to fill
- */
-void importNetgenParams(const std::string param_file, netgen_params& aParams){
+void importDefaultNetgenParams(const std::string param_file, netgen_params& aParams)
+{
   std::ifstream myfile(param_file);
   std::string line;
+  // set the default type!
+  aParams.myType = Hypo;
 
   std::getline(myfile, line);
   aParams.has_netgen_param = std::stoi(line);
@@ -131,7 +131,54 @@ void importNetgenParams(const std::string param_file, netgen_params& aParams){
   aParams.maxElementVolume = std::stod(line);
   std::getline(myfile, line);
   aParams.maxElementVolume = std::stoi(line);
+  myfile.close();
+}
+
+void importSimple2D3DNetgenParams(const std::string param_file, netgen_params& aParams, bool is3D )
+{
+  std::ifstream myfile(param_file);
+  std::string line;
+
+  aParams.myType = !is3D ? Simple2D : Simple3D;
+  std::getline(myfile, line);
+  aParams.has_netgen_param = std::stoi(line); // 1
+  std::getline(myfile, line);
+  aParams.numberOfSegments = std::stoi(line); // segments      (int)
+  std::getline(myfile, line);
+  aParams.localLength = std::stod(line);      // localLenght  (double)
+  std::getline(myfile, line);
+  aParams.maxElementArea = std::stod(line);   // max area     (double)
+  if ( is3D )
+  {
+    std::getline(myfile, line);
+    aParams.maxElementVol = std::stod(line); // max volume    (double)
+  }    
+  std::getline(myfile, line);
+  aParams.allowQuadrangles = std::stoi(line); // int
 
+  myfile.close();
+};
+
+/**
+ * @brief Import a param_file into a netgen_params structure
+ *
+ * @param param_file Name of the file
+ * @param aParams Structure to fill
+ */
+void importNetgenParams(const std::string param_file, netgen_params& aParams){
+  
+  if ( param_file.find("simple2D") != std::string::npos || param_file.find("simple3D") != std::string::npos /*support simple 2D && 3D*/ )
+  {
+    importSimple2D3DNetgenParams( param_file, aParams, bool(param_file.find("simple3D") != std::string::npos) );
+  }
+  else if ( param_file.find("maxarea") == std::string::npos && param_file.find("lenghtfromedge") == std::string::npos /*hypo file for 2D SA*/)
+  {
+    importDefaultNetgenParams( param_file, aParams );
+  }
+  else
+  {
+    aParams.has_netgen_param = false;
+  }
 };
 
 /**
@@ -141,32 +188,38 @@ void importNetgenParams(const std::string param_file, netgen_params& aParams){
  * @param aParams the object
  */
 void exportNetgenParams(const std::string param_file, netgen_params& aParams){
-  std::ofstream myfile(param_file);
-  myfile << aParams.has_netgen_param << std::endl;
-  myfile << aParams.maxh << std::endl;
-  myfile << aParams.minh << std::endl;
-  myfile << aParams.segmentsperedge << std::endl;
-  myfile << aParams.grading << std::endl;
-  myfile << aParams.curvaturesafety << std::endl;
-  myfile << aParams.secondorder << std::endl;
-  myfile << aParams.quad << std::endl;
-  myfile << aParams.optimize << std::endl;
-  myfile << aParams.fineness << std::endl;
-  myfile << aParams.uselocalh << std::endl;
-  myfile << aParams.merge_solids << std::endl;
-  myfile << aParams.chordalError << std::endl;
-  myfile << aParams.optsteps2d << std::endl;
-  myfile << aParams.optsteps3d << std::endl;
-  myfile << aParams.elsizeweight << std::endl;
-  myfile << aParams.opterrpow << std::endl;
-  myfile << aParams.delaunay << std::endl;
-  myfile << aParams.checkoverlap << std::endl;
-  myfile << aParams.checkchartboundary << std::endl;
-  myfile << aParams.closeedgefac << std::endl;
-  myfile << aParams.nbThreads << std::endl;
-  myfile << aParams.has_local_size << std::endl;
-  myfile << aParams.meshsizefilename << std::endl;
-  myfile << aParams.has_maxelementvolume_hyp << std::endl;
-  myfile << aParams.maxElementVolume << std::endl;
-  myfile << aParams.has_LengthFromEdges_hyp << std::endl;
+  if ( aParams.myType == Hypo ){
+    std::ofstream myfile(param_file);
+    myfile << aParams.has_netgen_param << std::endl;
+    myfile << aParams.maxh << std::endl;
+    myfile << aParams.minh << std::endl;
+    myfile << aParams.segmentsperedge << std::endl;
+    myfile << aParams.grading << std::endl;
+    myfile << aParams.curvaturesafety << std::endl;
+    myfile << aParams.secondorder << std::endl;
+    myfile << aParams.quad << std::endl;
+    myfile << aParams.optimize << std::endl;
+    myfile << aParams.fineness << std::endl;
+    myfile << aParams.uselocalh << std::endl;
+    myfile << aParams.merge_solids << std::endl;
+    myfile << aParams.chordalError << std::endl;
+    myfile << aParams.optsteps2d << std::endl;
+    myfile << aParams.optsteps3d << std::endl;
+    myfile << aParams.elsizeweight << std::endl;
+    myfile << aParams.opterrpow << std::endl;
+    myfile << aParams.delaunay << std::endl;
+    myfile << aParams.checkoverlap << std::endl;
+    myfile << aParams.checkchartboundary << std::endl;
+    myfile << aParams.closeedgefac << std::endl;
+    myfile << aParams.nbThreads << std::endl;
+    myfile << aParams.has_local_size << std::endl;
+    myfile << aParams.meshsizefilename << std::endl;
+    myfile << aParams.has_maxelementvolume_hyp << std::endl;
+    myfile << aParams.maxElementVolume << std::endl;
+    myfile << aParams.has_LengthFromEdges_hyp << std::endl;
+  }
+  else if ( aParams.myType == Simple2D )
+  {
+    // TODO: Export the 2D && 3D simple versions
+  }
 };
index 7fe0126373c1ef5e9c04eb27137a3bc6e67f777b..182c4ceaf91e72aaf52b2d0bb330a3164adac1df 100644 (file)
 
 #include <string>
 
+enum hypoType{
+  Hypo = 1, /*the default hypo common for 1D2D,1D2D3D && 3D*/
+  Simple2D,
+  Simple3D
+};
+
 struct netgen_params{
   // Params from NETGENPlugin_Mesher
   // True if _hypParameters is not null
   bool has_netgen_param=true;
+  hypoType myType;
   double maxh;
   double minh;
   double segmentsperedge;
@@ -70,11 +77,23 @@ struct netgen_params{
   // Params from NETGEN2D
   bool has_LengthFromEdges_hyp=false;
 
+  /////////////////////////////////////
+  // Quantities proper of Simple2D
+  int numberOfSegments; // maybe reuse segmentsperedge (?)
+  double localLength;
+  double maxElementArea;
+  bool allowQuadrangles;
+  //// Quantities proper of Simple3D
+  double maxElementVol;
 };
 
 void printNetgenParams(netgen_params& aParams);
 
 void importNetgenParams(const std::string param_file, netgen_params& aParams);
+void importDefaultNetgenParams(const std::string param_file, netgen_params& aParams);
+void importSimple2D3DNetgenParams(const std::string param_file, netgen_params& aParams, bool is3D );
 void exportNetgenParams(const std::string param_file, netgen_params& aParams);
+// TODO symple param to be used with netgen1d2d
+
 
 #endif
index 34cf316240adaa34f040405fa7a8c43d07b48cbe..9f7268332c7f910e65dd2b45167483e8d3567121 100644 (file)
@@ -2270,7 +2270,6 @@ NETGENPlugin_Mesher::AddSegmentsToMesh(netgen::Mesh&                    ngMesh,
 
   const bool wasNgMeshEmpty = ( ngMesh.GetNP() < 1 ); /* true => this method is called by
                                                          NETGENPlugin_NETGEN_2D_ONLY */
-
   // map for nodes on vertices since they can be shared between wires
   // ( issue 0020676, face_int_box.brep) and nodes built by NETGEN
   map<const SMDS_MeshNode*, int > node2ngID;
@@ -2330,8 +2329,10 @@ NETGENPlugin_Mesher::AddSegmentsToMesh(netgen::Mesh&                    ngMesh,
         continue;
 
       int ngID1 = ngMesh.GetNP() + 1, ngID2 = ngID1+1;
+
       if ( onVertex || ( !wasNgMeshEmpty && onEdge ) || helper.IsRealSeam( posShapeID ))
         ngID1 = node2ngID.insert( make_pair( n, ngID1 )).first->second;
+      
       if ( ngID1 > ngMesh.GetNP() )
       {
         netgen::MeshPoint mp( netgen::Point<3> (n->X(), n->Y(), n->Z()) );
@@ -2360,11 +2361,9 @@ NETGENPlugin_Mesher::AddSegmentsToMesh(netgen::Mesh&                    ngMesh,
       for ( int iEnd = 0; iEnd < 2; ++iEnd)
       {
         const UVPtStruct& pnt = uvPtVec[ i + iEnd ];
-
         seg.epgeominfo[ iEnd ].dist = pnt.param; // param on curve
         seg.epgeominfo[ iEnd ].u    = pnt.u;
         seg.epgeominfo[ iEnd ].v    = pnt.v;
-
         // find out edge id and node parameter on edge
         onVertex = ( pnt.normParam + 1e-10 > vertexNormPar );
         if ( onVertex || posShapeID != posID )
@@ -2495,7 +2494,6 @@ int NETGENPlugin_Mesher::FillSMesh(const netgen::OCCGeometry&          occgeo,
   int nbSeg = ngMesh.GetNSeg();
   int nbFac = ngMesh.GetNSE();
   int nbVol = ngMesh.GetNE();
-
   SMESHDS_Mesh* meshDS = sMesh.GetMeshDS();
 
   // quadHelper is used for either
@@ -2879,528 +2877,787 @@ namespace
   const double volOptimizeTime = 0.77;
 }
 
-//=============================================================================
-/*!
- * Here we are going to use the NETGEN mesher
- */
-//=============================================================================
-
-bool NETGENPlugin_Mesher::Compute()
+int NETGENPlugin_Mesher::FillInternalElements( NETGENPlugin_NetgenLibWrapper& ngLib, NETGENPlugin_Internals& internals, netgen::OCCGeometry& occgeo, 
+                                                NETGENPlugin_ngMeshInfo& initState, SMESH_MesherHelper &quadHelper, list< SMESH_subMesh* >* meshedSM )
 {
-  NETGENPlugin_NetgenLibWrapper ngLib;
+  SMESH_Comment comment;
+  int err;
+  int startWith = netgen::MESHCONST_ANALYSE;
+  int endWith   = netgen::MESHCONST_ANALYSE;
+
+  // load internal shapes into OCCGeometry
+  netgen::OCCGeometry intOccgeo;
+  internals.getInternalEdges( intOccgeo.fmap, intOccgeo.emap, intOccgeo.vmap, meshedSM );
+  intOccgeo.boundingbox = occgeo.boundingbox;
+  intOccgeo.shape = occgeo.shape;
+  intOccgeo.face_maxh.SetSize(intOccgeo.fmap.Extent());
+  intOccgeo.face_maxh = netgen::mparam.maxh;
+  netgen::Mesh *tmpNgMesh = NULL;  
+
+  try
+  {
+    OCC_CATCH_SIGNALS;
+    // compute local H on internal shapes in the main mesh
+    //OCCSetLocalMeshSize(intOccgeo, *_ngMesh); it deletes _ngMesh->localH
 
-  netgen::MeshingParameters& mparams = netgen::mparam;
+    // let netgen create a temporary mesh
+    ngLib.GenerateMesh(intOccgeo, startWith, endWith, tmpNgMesh);
 
-  SMESH_ComputeErrorPtr error = SMESH_ComputeError::New();
-  SMESH_MesherHelper quadHelper( *_mesh );
-  quadHelper.SetIsQuadratic( mparams.secondorder );
+    if ( netgen::multithread.terminate )
+      return false;
 
-  // -------------------------
-  // Prepare OCC geometry
-  // -------------------------
+    // copy LocalH from the main to temporary mesh
+    initState.transferLocalH( _ngMesh, tmpNgMesh );
 
-  netgen::OCCGeometry occgeo;
-  list< SMESH_subMesh* > meshedSM[3]; // for 0-2 dimensions
-  NETGENPlugin_Internals internals( *_mesh, _shape, _isVolume );
-  PrepareOCCgeometry( occgeo, _shape, *_mesh, meshedSM, &internals );
-  _occgeom = &occgeo;
+    // compute mesh on internal edges
+    startWith = endWith = netgen::MESHCONST_MESHEDGES;
+    err = ngLib.GenerateMesh(intOccgeo, startWith, endWith, tmpNgMesh);
 
-  _totalTime = edgeFaceMeshingTime;
-  if ( _optimize )
-    _totalTime += faceOptimizTime;
-  if ( _isVolume )
-    _totalTime += voluMeshingTime + ( _optimize ? volOptimizeTime : 0 );
-  double doneTime = 0;
-  _ticTime = -1;
-  _progressTic = 1;
-  _curShapeIndex = -1;
+    comment << text(err);
+  }
+  catch (Standard_Failure& ex)
+  {
+    comment << text(ex);
+    err = 1;
+  }
+  initState.restoreLocalH( tmpNgMesh );
 
-  // -------------------------
-  // Generate the mesh
-  // -------------------------
+  // fill SMESH by netgen mesh
+  vector< const SMDS_MeshNode* > tmpNodeVec;
+  FillSMesh( intOccgeo, *tmpNgMesh, initState, *_mesh, tmpNodeVec, comment, &quadHelper );
+  err = ( err || !comment.empty() );
 
-  _ngMesh = NULL;
-  NETGENPlugin_ngMeshInfo initState; // it remembers size of ng mesh equal to size of Smesh
+  nglib::Ng_DeleteMesh((nglib::Ng_Mesh*)tmpNgMesh);
+  return err;
+}
 
+bool NETGENPlugin_Mesher::Fill2DViscousLayer( netgen::OCCGeometry& occgeo, vector< const SMDS_MeshNode* >& nodeVec, 
+                                              NETGENPlugin_Internals* internals, NETGENPlugin_ngMeshInfo& initState )
+{
+  SMESH_ComputeErrorPtr error = SMESH_ComputeError::New();
   SMESH_Comment comment;
-  int err = 0;
-
-  // vector of nodes in which node index == netgen ID
-  vector< const SMDS_MeshNode* > nodeVec;
 
+  bool toCompute = _isViscousLayers2D ||
+                    ( !occgeo.fmap.IsEmpty() && 
+                        StdMeshers_ViscousLayers2D::HasProxyMesh( TopoDS::Face( occgeo.fmap(1) ), *_mesh )
+                    );
+  if ( toCompute )
   {
-    // ----------------
-    // compute 1D mesh
-    // ----------------
-    if ( _simpleHyp )
+    if ( !internals->hasInternalVertexInFace() ) {
+      FillSMesh( occgeo, *_ngMesh, initState, *_mesh, nodeVec, comment );
+      initState = NETGENPlugin_ngMeshInfo(_ngMesh);
+    }
+    SMESH_ProxyMesh::Ptr viscousMesh;
+    SMESH_MesherHelper   helper( *_mesh );
+    for ( int faceID = 1; faceID <= occgeo.fmap.Extent(); ++faceID )
     {
-      // not to RestrictLocalH() according to curvature during MESHCONST_ANALYSE
-      mparams.uselocalh = false;
-      mparams.grading = 0.8; // not limitited size growth
+      const TopoDS_Face& F = TopoDS::Face( occgeo.fmap( faceID ));
+      viscousMesh = StdMeshers_ViscousLayers2D::Compute( *_mesh, F );
+      if ( !viscousMesh )
+        return false;
+      if ( viscousMesh->NbProxySubMeshes() == 0 )
+        continue;
+      // exclude from computation ng segments built on EDGEs of F
+      for (int i = 1; i <= _ngMesh->GetNSeg(); i++)
+      {
+        netgen::Segment & seg = _ngMesh->LineSegment(i);
+        if (seg.si == faceID)
+          seg.si = 0;
+      }
+      // add new segments to _ngMesh instead of excluded ones
+      helper.SetSubShape( F );
+      TSideVector wires = StdMeshers_FaceSide::GetFaceWires( F, *_mesh, /*skipMediumNodes=*/true, error, &helper, viscousMesh );
+      error = AddSegmentsToMesh( *_ngMesh, occgeo, wires, helper, nodeVec );
 
-      if ( _simpleHyp->GetNumberOfSegments() )
-        // nb of segments
-        mparams.maxh = occgeo.boundingbox.Diam();
-      else
-        // segment length
-        mparams.maxh = _simpleHyp->GetLocalLength();
+      if ( !error ) error = SMESH_ComputeError::New();
     }
+    initState = NETGENPlugin_ngMeshInfo(_ngMesh);
+  }
+  return true;
+}
 
-    if ( mparams.maxh == 0.0 )
-      mparams.maxh = occgeo.boundingbox.Diam();
-    if ( _simpleHyp || ( mparams.minh == 0.0 && _fineness != NETGENPlugin_Hypothesis::UserDefined))
-      mparams.minh = GetDefaultMinSize( _shape, mparams.maxh );
+bool NETGENPlugin_Mesher::Fill3DViscousLayerAndQuadAdaptor( netgen::OCCGeometry& occgeo, vector< const SMDS_MeshNode* >& nodeVec, 
+                                                            netgen::MeshingParameters &mparams, NETGENPlugin_ngMeshInfo& initState,
+                                                            list< SMESH_subMesh* >* meshedSM, SMESH_MesherHelper &quadHelper, int& err )
+{
+  SMESH_Comment comment;
+  if ( _isVolume && ( !meshedSM[ MeshDim_2D ].empty() || mparams.quad || _viscousLayersHyp ) )
+  {
+    // load SMESH with computed segments and faces
+    FillSMesh( occgeo, *_ngMesh, initState, *_mesh, nodeVec, comment, &quadHelper );
+
+    // compute prismatic boundary volumes
+    smIdType nbQuad = _mesh->NbQuadrangles();
+    SMESH_ProxyMesh::Ptr viscousMesh;
+    if ( _viscousLayersHyp )
+    {
+      viscousMesh = _viscousLayersHyp->Compute( *_mesh, _shape );
+      if ( !viscousMesh )
+        return false;
+    }
+    // compute pyramids on quadrangles
+    vector<SMESH_ProxyMesh::Ptr> pyramidMeshes( occgeo.somap.Extent() );
+    if ( nbQuad > 0 )
+      for ( int iS = 1; iS <= occgeo.somap.Extent(); ++iS )
+      {
+        StdMeshers_QuadToTriaAdaptor* adaptor = new StdMeshers_QuadToTriaAdaptor;
+        pyramidMeshes[ iS-1 ].reset( adaptor );
+        bool ok = adaptor->Compute( *_mesh, occgeo.somap(iS), viscousMesh.get() );
+        if ( !ok )
+          return false;
+      }
+    
+    // add proxy faces to NG mesh
+    list< SMESH_subMesh* > viscousSM;
+    for ( int iS = 1; iS <= occgeo.somap.Extent(); ++iS )
+    {
+      list< SMESH_subMesh* > quadFaceSM;
+      for (TopExp_Explorer face(occgeo.somap(iS), TopAbs_FACE); face.More(); face.Next())
+        if ( pyramidMeshes[iS-1] && pyramidMeshes[iS-1]->GetProxySubMesh( face.Current() ))
+        {
+          quadFaceSM.push_back( _mesh->GetSubMesh( face.Current() ));
+          meshedSM[ MeshDim_2D ].remove( quadFaceSM.back() );
+        }
+        else if ( viscousMesh && viscousMesh->GetProxySubMesh( face.Current() ))
+        {
+          viscousSM.push_back( _mesh->GetSubMesh( face.Current() ));
+          meshedSM[ MeshDim_2D ].remove( viscousSM.back() );
+        }
+      if ( !quadFaceSM.empty() )
+        FillNgMesh(occgeo, *_ngMesh, nodeVec, quadFaceSM, &quadHelper, pyramidMeshes[iS-1]);
+    }
+    if ( !viscousSM.empty() )
+      FillNgMesh(occgeo, *_ngMesh, nodeVec, viscousSM, &quadHelper, viscousMesh );
+
+    // fill _ngMesh with faces of sub-meshes
+    err = ! ( FillNgMesh(occgeo, *_ngMesh, nodeVec, meshedSM[ MeshDim_2D ], &quadHelper));
+    initState = NETGENPlugin_ngMeshInfo(_ngMesh, /*checkRemovedElems=*/true);
+  }
+  return true;
+}
+
+void NETGENPlugin_Mesher::CallNetgenConstAnalysis( NETGENPlugin_NetgenLibWrapper& ngLib, netgen::MeshingParameters &mparams, netgen::OCCGeometry& occgeo )
+{
+  SMESH_Comment comment;
+  int err;  
+  int startWith = netgen::MESHCONST_ANALYSE;
+  int endWith   = netgen::MESHCONST_ANALYSE;
+  try
+  {
+    OCC_CATCH_SIGNALS;
+
+    err = ngLib.GenerateMesh(occgeo, startWith, endWith, _ngMesh );
+    // if(netgen::multithread.terminate)
+    //   return false;
+    comment << text(err);
+  }
+  catch (Standard_Failure& ex)
+  {
+    comment << text(ex);
+  }
+  catch (netgen::NgException & ex)
+  {
+    comment << text(ex);
+#ifdef NETGEN_V6
+    bool hasSizeFile = !mparams.meshsizefilename.empty();
+#else
+    bool hasSizeFile = mparams.meshsizefilename;
+#endif
+    if ( hasSizeFile )
+      throw SMESH_ComputeError(COMPERR_BAD_PARMETERS, comment );
+  }
+
+  ngLib.setMesh(( Ng_Mesh*) _ngMesh );
+}
+
+int NETGENPlugin_Mesher::CallNetgenMeshEdges( NETGENPlugin_NetgenLibWrapper& ngLib, netgen::OCCGeometry& occgeo )
+{
+  SMESH_Comment comment;
+  int err = 0;  
+  int startWith = netgen::MESHCONST_MESHEDGES; 
+  int endWith   = netgen::MESHCONST_MESHEDGES;
+  try
+  {
+    OCC_CATCH_SIGNALS;
+
+    err = ngLib.GenerateMesh(occgeo, startWith, endWith, _ngMesh );
+    // if(netgen::multithread.terminate)
+    //   return false;
+    comment << text(err);
+  }
+  catch (Standard_Failure& ex)
+  {
+    comment << text(ex);
+    err = 1;
+  }
+  catch (netgen::NgException & ex)
+  {
+    comment << text(ex);
+    err = 1;
+  }
+  return err;
+}
+
+int NETGENPlugin_Mesher::CallNetgenMeshFaces( NETGENPlugin_NetgenLibWrapper& ngLib, netgen::OCCGeometry& occgeo, SMESH_Comment& comment )
+{
+  int err = 0;  
+  int startWith = netgen::MESHCONST_MESHSURFACE; 
+  int endWith   =  _optimize ? netgen::MESHCONST_OPTSURFACE : netgen::MESHCONST_MESHSURFACE;
+  try
+  {
+    OCC_CATCH_SIGNALS;
+
+    err = ngLib.GenerateMesh(occgeo, startWith, endWith, _ngMesh );
+    // if(netgen::multithread.terminate)
+    //   return false;
+    comment << text(err);
+  }
+  catch (Standard_Failure& ex)
+  {
+    comment << text(ex);
+  }
+  catch (netgen::NgException & ex)
+  {
+    comment << text(ex);
+  }
+  return err;
+}
+
+int NETGENPlugin_Mesher::CallNetgenMeshVolumens( NETGENPlugin_NetgenLibWrapper& ngLib, netgen::OCCGeometry& occgeo, SMESH_Comment& comment )
+{
+  // Let netgen compute 3D mesh
+  int err = 0;
+  int startWith = netgen::MESHCONST_MESHVOLUME;
+  int endWith   = netgen::MESHCONST_MESHVOLUME;
+  try
+  {
+    OCC_CATCH_SIGNALS;
+
+    err = ngLib.GenerateMesh(occgeo, startWith, endWith);
 
-    // Local size on faces
-    occgeo.face_maxh = mparams.maxh;
+    if ( netgen::multithread.terminate )
+      return false;
+
+    if ( comment.empty() ) // do not overwrite a previous error
+      comment << text(err);
+  }
+  catch (Standard_Failure& ex)
+  {
+    if ( comment.empty() ) // do not overwrite a previous error
+      comment << text(ex);
+    err = 1;
+  }
+  catch (netgen::NgException& exc)
+  {
+    if ( comment.empty() ) // do not overwrite a previous error
+      comment << text(exc);
+    err = 1;
+  }
+  // _ticTime = ( doneTime += voluMeshingTime ) / _totalTime / _progressTic;
 
-    // Let netgen create _ngMesh and calculate element size on not meshed shapes
-    int startWith = netgen::MESHCONST_ANALYSE;
-    int endWith   = netgen::MESHCONST_ANALYSE;
+  // Let netgen optimize 3D mesh
+  if ( !err && _optimize )
+  {
+    startWith = endWith = netgen::MESHCONST_OPTVOLUME;
     try
     {
       OCC_CATCH_SIGNALS;
 
-      err = ngLib.GenerateMesh(occgeo, startWith, endWith, _ngMesh);
+      err = ngLib.GenerateMesh(occgeo, startWith, endWith);
 
-      if(netgen::multithread.terminate)
+      if ( netgen::multithread.terminate )
         return false;
 
-      comment << text(err);
+      if ( comment.empty() ) // do not overwrite a previous error
+        comment << text(err);
     }
     catch (Standard_Failure& ex)
     {
-      comment << text(ex);
+      if ( comment.empty() ) // do not overwrite a previous error
+        comment << text(ex);
     }
-    catch (netgen::NgException & ex)
+    catch (netgen::NgException& exc)
     {
-      comment << text(ex);
-#ifdef NETGEN_V6
-      bool hasSizeFile = !mparams.meshsizefilename.empty();
-#else
-      bool hasSizeFile = mparams.meshsizefilename;
-#endif
-      if ( hasSizeFile )
-        throw SMESH_ComputeError(COMPERR_BAD_PARMETERS, comment );
+      if ( comment.empty() ) // do not overwrite a previous error
+        comment << text(exc);
     }
-    err = 0; //- MESHCONST_ANALYSE isn't so important step
-    if ( !_ngMesh )
-      return false;
-    ngLib.setMesh(( Ng_Mesh*) _ngMesh );
-
-    _ngMesh->ClearFaceDescriptors(); // we make descriptors our-self
-
-    if ( !mparams.uselocalh ) // mparams.grading is not taken into account yet
-      _ngMesh->LocalHFunction().SetGrading( mparams.grading );
+  }
+  return err;
+}
 
-    if ( _simpleHyp )
+void NETGENPlugin_Mesher::MakeSecondOrder( NETGENPlugin_NetgenLibWrapper& ngLib, netgen::MeshingParameters &mparams, netgen::OCCGeometry& occgeo, 
+                                           list< SMESH_subMesh* >* meshedSM, NETGENPlugin_ngMeshInfo& initState, SMESH_Comment& comment )
+{
+  if ( mparams.secondorder > 0 )
+  {
+    try
     {
-      // Pass 1D simple parameters to NETGEN
-      // --------------------------------
-      double nbSeg   = (double) _simpleHyp->GetNumberOfSegments();
-      double segSize = _simpleHyp->GetLocalLength();
-      for ( int iE = 1; iE <= occgeo.emap.Extent(); ++iE )
+      OCC_CATCH_SIGNALS;
+      if ( !meshedSM[ MeshDim_1D ].empty() )
       {
-        const TopoDS_Edge& e = TopoDS::Edge( occgeo.emap(iE));
-        if ( nbSeg )
-          segSize = SMESH_Algo::EdgeLength( e ) / ( nbSeg - 0.4 );
-        setLocalSize( e, segSize, *_ngMesh );
+        // remove segments not attached to geometry (IPAL0052479)
+        for (int i = 1; i <= _ngMesh->GetNSeg(); ++i)
+        {
+          const netgen::Segment & seg = _ngMesh->LineSegment (i);
+          if ( seg.epgeominfo[ 0 ].edgenr == 0 )
+          {
+            _ngMesh->DeleteSegment( i );
+            initState._nbSegments--;
+          }
+        }
+        _ngMesh->Compress();
       }
+      // convert to quadratic
+#ifdef NETGEN_V6
+      occgeo.GetRefinement().MakeSecondOrder(*_ngMesh);
+#else
+      netgen::OCCRefinementSurfaces(occgeo).MakeSecondOrder(*_ngMesh);
+#endif
+
+      // care of elements already loaded to SMESH
+      // if ( initState._nbSegments > 0 )
+      //   makeQuadratic( occgeo.emap, _mesh );
+      // if ( initState._nbFaces > 0 )
+      //   makeQuadratic( occgeo.fmap, _mesh );
     }
-    else // if ( ! _simpleHyp )
+    catch (Standard_Failure& ex)
     {
-      // Local size on shapes
-      SetLocalSize( occgeo, *_ngMesh );
-      SetLocalSizeForChordalError( occgeo, *_ngMesh );
+      if ( comment.empty() ) // do not overwrite a previous error
+        comment << "Exception in netgen at passing to 2nd order ";
     }
-
-    // Precompute internal edges (issue 0020676) in order to
-    // add mesh on them correctly (twice) to netgen mesh
-    if ( !err && internals.hasInternalEdges() )
+    catch (netgen::NgException& exc)
     {
-      // load internal shapes into OCCGeometry
-      netgen::OCCGeometry intOccgeo;
-      internals.getInternalEdges( intOccgeo.fmap, intOccgeo.emap, intOccgeo.vmap, meshedSM );
-      intOccgeo.boundingbox = occgeo.boundingbox;
-      intOccgeo.shape = occgeo.shape;
-      intOccgeo.face_maxh.SetSize(intOccgeo.fmap.Extent());
-      intOccgeo.face_maxh = netgen::mparam.maxh;
-      netgen::Mesh *tmpNgMesh = NULL;
-      try
-      {
-        OCC_CATCH_SIGNALS;
-        // compute local H on internal shapes in the main mesh
-        //OCCSetLocalMeshSize(intOccgeo, *_ngMesh); it deletes _ngMesh->localH
+      if ( comment.empty() ) // do not overwrite a previous error
+        comment << exc.What();
+    }    
+  }    
+}
 
-        // let netgen create a temporary mesh
-        ngLib.GenerateMesh(intOccgeo, startWith, endWith, tmpNgMesh);
+// Define minimal parameter for netgen mesher
+void NETGENPlugin_Mesher::SetBasicMeshParameters( NETGENPlugin_NetgenLibWrapper& ngLib, netgen::MeshingParameters &mparams, netgen::OCCGeometry& occgeo )
+{
+  if ( _simpleHyp )
+  {
+    // not to RestrictLocalH() according to curvature during MESHCONST_ANALYSE
+    mparams.uselocalh = false;
+    mparams.grading = 0.8; // not limitited size growth
 
-        if ( netgen::multithread.terminate )
-          return false;
+    if ( _simpleHyp->GetNumberOfSegments() )
+      // nb of segments
+      mparams.maxh = occgeo.boundingbox.Diam();
+    else
+      // segment length
+      mparams.maxh = _simpleHyp->GetLocalLength();
+  }
 
-        // copy LocalH from the main to temporary mesh
-        initState.transferLocalH( _ngMesh, tmpNgMesh );
+  if ( mparams.maxh == 0.0 )
+    mparams.maxh = occgeo.boundingbox.Diam();
+  if ( _simpleHyp || ( mparams.minh == 0.0 && _fineness != NETGENPlugin_Hypothesis::UserDefined))
+    mparams.minh = GetDefaultMinSize( _shape, mparams.maxh );
 
-        // compute mesh on internal edges
-        startWith = endWith = netgen::MESHCONST_MESHEDGES;
-        err = ngLib.GenerateMesh(intOccgeo, startWith, endWith, tmpNgMesh);
+  // Local size on faces
+  occgeo.face_maxh = mparams.maxh;
 
-        comment << text(err);
-      }
-      catch (Standard_Failure& ex)
-      {
-        comment << text(ex);
-        err = 1;
-      }
-      initState.restoreLocalH( tmpNgMesh );
+  CallNetgenConstAnalysis( ngLib, mparams, occgeo ); // -> Initialize _ngMesh and set ngLib->_ngMesh = _ngMesh
 
-      // fill SMESH by netgen mesh
-      vector< const SMDS_MeshNode* > tmpNodeVec;
-      FillSMesh( intOccgeo, *tmpNgMesh, initState, *_mesh, tmpNodeVec, comment, &quadHelper );
-      err = ( err || !comment.empty() );
+  _ngMesh->ClearFaceDescriptors(); // we make descriptors our-self
 
-      nglib::Ng_DeleteMesh((nglib::Ng_Mesh*)tmpNgMesh);
-    }
+  if ( !mparams.uselocalh ) // mparams.grading is not taken into account yet
+    _ngMesh->LocalHFunction().SetGrading( mparams.grading );
 
-    // Fill _ngMesh with nodes and segments of computed submeshes
-    if ( !err )
+  if ( _simpleHyp )
+  {
+    // Pass 1D simple parameters to NETGEN
+    // --------------------------------
+    double nbSeg   = (double) _simpleHyp->GetNumberOfSegments();
+    double segSize = _simpleHyp->GetLocalLength();
+    for ( int iE = 1; iE <= occgeo.emap.Extent(); ++iE )
     {
-      err = ! ( FillNgMesh(occgeo, *_ngMesh, nodeVec, meshedSM[ MeshDim_0D ]) &&
-                FillNgMesh(occgeo, *_ngMesh, nodeVec, meshedSM[ MeshDim_1D ], &quadHelper));
+      const TopoDS_Edge& e = TopoDS::Edge( occgeo.emap(iE));
+      if ( nbSeg )
+        segSize = SMESH_Algo::EdgeLength( e ) / ( nbSeg - 0.4 );
+      setLocalSize( e, segSize, *_ngMesh );
     }
-    initState = NETGENPlugin_ngMeshInfo(_ngMesh);
-
-    // Compute 1d mesh
-    if (!err)
-    {
-      startWith = endWith = netgen::MESHCONST_MESHEDGES;
-      try
-      {
-        OCC_CATCH_SIGNALS;
-
-        err = ngLib.GenerateMesh(occgeo, startWith, endWith);
+  }
+  else // if ( ! _simpleHyp )
+  {
+    // Local size on shapes
+    SetLocalSize( occgeo, *_ngMesh );
+    SetLocalSizeForChordalError( occgeo, *_ngMesh );
+  }
+}
 
-        if ( netgen::multithread.terminate )
-          return false;
+void NETGENPlugin_Mesher::SetBasicMeshParametersFor2D( netgen::OCCGeometry& occgeo, vector< const SMDS_MeshNode* >& nodeVec, 
+                                                        netgen::MeshingParameters &mparams, NETGENPlugin_Internals* internals, 
+                                                          NETGENPlugin_ngMeshInfo& initState )
+{
+  SMESH_Comment comment;
+  mparams.uselocalh = true;
 
-        comment << text(err);
+  if ( _simpleHyp ) {
+    if ( double area = _simpleHyp->GetMaxElementArea() ) {
+      // face area
+      mparams.maxh = sqrt(2. * area/sqrt(3.0));
+      mparams.grading = 0.4; // moderate size growth
+    }
+    else {
+      // length from edges
+      if ( _ngMesh->GetNSeg() ) {
+        double edgeLength = 0;
+        TopTools_MapOfShape visitedEdges;
+        for ( TopExp_Explorer exp( _shape, TopAbs_EDGE ); exp.More(); exp.Next() )
+          if( visitedEdges.Add(exp.Current()) )
+            edgeLength += SMESH_Algo::EdgeLength( TopoDS::Edge( exp.Current() ));
+        // we have to multiply length by 2 since for each TopoDS_Edge there
+        // are double set of NETGEN edges, in other words, we have to
+        // divide _ngMesh->GetNSeg() by 2.
+        mparams.maxh = 2*edgeLength / _ngMesh->GetNSeg();
       }
-      catch (Standard_Failure& ex)
-      {
-        comment << text(ex);
-        err = 1;
+      else {
+        mparams.maxh = 1000;
       }
+      mparams.grading = 0.2; // slow size growth
     }
-    if ( _isVolume )
-      _ticTime = ( doneTime += edgeMeshingTime ) / _totalTime / _progressTic;
+    mparams.quad = _simpleHyp->GetAllowQuadrangles();
+    mparams.maxh = min( mparams.maxh, occgeo.boundingbox.Diam()/2 );
+    _ngMesh->SetGlobalH (mparams.maxh);
+    netgen::Box<3> bb = occgeo.GetBoundingBox();
+    bb.Increase (bb.Diam()/20);
+    _ngMesh->SetLocalH (bb.PMin(), bb.PMax(), mparams.grading);
+  }
 
-    mparams.uselocalh = true; // restore as it is used at surface optimization
+  // Care of vertices internal in faces (issue 0020676)
+  if ( internals->hasInternalVertexInFace() )
+  {
+    // store computed segments in SMESH in order not to create SMESH
+    // edges for ng segments added by AddIntVerticesInFaces()
+    FillSMesh( occgeo, *_ngMesh, initState, *_mesh, nodeVec, comment );
+    // add segments to faces with internal vertices
+    AddIntVerticesInFaces( occgeo, *_ngMesh, nodeVec, *internals );
+    initState = NETGENPlugin_ngMeshInfo(_ngMesh);
+  }
+}
 
-    // ---------------------
-    // compute surface mesh
-    // ---------------------
-    if (!err)
-    {
-      // Pass 2D simple parameters to NETGEN
-      if ( _simpleHyp ) {
-        if ( double area = _simpleHyp->GetMaxElementArea() ) {
-          // face area
-          mparams.maxh = sqrt(2. * area/sqrt(3.0));
-          mparams.grading = 0.4; // moderate size growth
-        }
-        else {
-          // length from edges
-          if ( _ngMesh->GetNSeg() ) {
-            double edgeLength = 0;
-            TopTools_MapOfShape visitedEdges;
-            for ( TopExp_Explorer exp( _shape, TopAbs_EDGE ); exp.More(); exp.Next() )
-              if( visitedEdges.Add(exp.Current()) )
-                edgeLength += SMESH_Algo::EdgeLength( TopoDS::Edge( exp.Current() ));
-            // we have to multiply length by 2 since for each TopoDS_Edge there
-            // are double set of NETGEN edges, in other words, we have to
-            // divide _ngMesh->GetNSeg() by 2.
-            mparams.maxh = 2*edgeLength / _ngMesh->GetNSeg();
-          }
-          else {
-            mparams.maxh = 1000;
-          }
-          mparams.grading = 0.2; // slow size growth
-        }
-        mparams.quad = _simpleHyp->GetAllowQuadrangles();
-        mparams.maxh = min( mparams.maxh, occgeo.boundingbox.Diam()/2 );
-        _ngMesh->SetGlobalH (mparams.maxh);
-        netgen::Box<3> bb = occgeo.GetBoundingBox();
-        bb.Increase (bb.Diam()/20);
-        _ngMesh->SetLocalH (bb.PMin(), bb.PMax(), mparams.grading);
-      }
+void NETGENPlugin_Mesher::SetBasicMeshParametersFor3D( NETGENPlugin_NetgenLibWrapper& ngLib, netgen::OCCGeometry& occgeo, 
+                                                        vector< const SMDS_MeshNode* >& nodeVec, netgen::MeshingParameters &mparams, 
+                                                          NETGENPlugin_Internals* internals, NETGENPlugin_ngMeshInfo& initState,  SMESH_MesherHelper &quadHelper, 
+                                                            SMESH_Comment& comment )
+{
+  const NETGENPlugin_SimpleHypothesis_3D* simple3d = dynamic_cast< const NETGENPlugin_SimpleHypothesis_3D* > ( _simpleHyp );
+  if ( simple3d ) {
+    _ngMesh->Compress();
+    if ( double vol = simple3d->GetMaxElementVolume() ) {
+      // max volume
+      mparams.maxh = pow( 72, 1/6. ) * pow( vol, 1/3. );
+      mparams.maxh = min( mparams.maxh, occgeo.boundingbox.Diam()/2 );
+    }
+    else {
+      // length from faces
+      mparams.maxh = _ngMesh->AverageH();
+    }
+    _ngMesh->SetGlobalH (mparams.maxh);
+    mparams.grading = 0.4;
+    ngLib.CalcLocalH( ngLib._ngMesh );
+  }
+  // Care of vertices internal in solids and internal faces (issue 0020676)
+  if ( internals->hasInternalVertexInSolid() || internals->hasInternalFaces() )
+  {
+    // store computed faces in SMESH in order not to create SMESH
+    // faces for ng faces added here
+    FillSMesh( occgeo, *_ngMesh, initState, *_mesh, nodeVec, comment, &quadHelper );
+    // add ng faces to solids with internal vertices
+    AddIntVerticesInSolids( occgeo, *_ngMesh, nodeVec, *internals );
+    // duplicate mesh faces on internal faces
+    FixIntFaces( occgeo, *_ngMesh, *internals );
+    initState = NETGENPlugin_ngMeshInfo(_ngMesh);
+  }
+}
 
-      // Care of vertices internal in faces (issue 0020676)
-      if ( internals.hasInternalVertexInFace() )
-      {
-        // store computed segments in SMESH in order not to create SMESH
-        // edges for ng segments added by AddIntVerticesInFaces()
-        FillSMesh( occgeo, *_ngMesh, initState, *_mesh, nodeVec, comment );
-        // add segments to faces with internal vertices
-        AddIntVerticesInFaces( occgeo, *_ngMesh, nodeVec, internals );
-        initState = NETGENPlugin_ngMeshInfo(_ngMesh);
-      }
+int NETGENPlugin_Mesher::Fill0D1DElements( netgen::OCCGeometry& occgeo, vector< const SMDS_MeshNode* >& nodeVec, list< SMESH_subMesh* >* meshedSM, SMESH_MesherHelper &quadHelper )
+{
+  int err;
+  err = ! ( FillNgMesh( occgeo, *_ngMesh, nodeVec, meshedSM[ MeshDim_0D ]) &&
+            FillNgMesh( occgeo, *_ngMesh, nodeVec, meshedSM[ MeshDim_1D ], &quadHelper));
+  return err;
+}
+
+void NETGENPlugin_Mesher::InitialSetup( NETGENPlugin_NetgenLibWrapper& ngLib, netgen::OCCGeometry& occgeo, 
+                                        list< SMESH_subMesh* >* meshedSM, NETGENPlugin_Internals* internals,
+                                        SMESH_MesherHelper &quadHelper, NETGENPlugin_ngMeshInfo& initState, 
+                                        netgen::MeshingParameters &mparams )
+{
+  int err = 0;
+  // Init occ geometry maps for non meshed object and fill meshedSM with premeshed objects
+  PrepareOCCgeometry( occgeo, _shape, *_mesh, meshedSM, internals );
+  _occgeom = &occgeo;
+  _ngMesh = NULL;
 
-      // Build viscous layers
-      if (( _isViscousLayers2D ) ||
-          ( !occgeo.fmap.IsEmpty() &&
-            StdMeshers_ViscousLayers2D::HasProxyMesh( TopoDS::Face( occgeo.fmap(1) ), *_mesh )))
+  // Fill basic mesh param and 
+  // define _ngMesh
+  SetBasicMeshParameters( ngLib, mparams, occgeo );
+  // Mesh internal faces and edges with netgen and then fill smesh with those entries!
+  if ( internals->hasInternalEdges() )
+    FillInternalElements( ngLib, *internals, occgeo, initState, quadHelper, meshedSM );
+   
+}
+
+void NETGENPlugin_Mesher::InitialSetupSA( NETGENPlugin_NetgenLibWrapper& ngLib, netgen::OCCGeometry& occgeo, 
+                                          list< SMESH_subMesh* >* meshedSM, NETGENPlugin_Internals* internals,
+                                          SMESH_MesherHelper &quadHelper, NETGENPlugin_ngMeshInfo& initState, 
+                                          netgen::MeshingParameters &mparams, bool useFMapFunction )
+{
+  // Check if *_mesh support submeshing so the geometry associated to it is used 
+  // so the original version of PrepareOCCgeometry is used.
+  // bool getSubmesh = false;
+  // for ( TopoDS_Iterator it( _mesh->GetShapeToMesh() ); it.More(); it.Next() )
+  // {
+  //   getSubmesh = _mesh->GetSubMesh( it.Value() );
+  //   if ( getSubmesh ) break;
+  // }
+  
+  // if ( !getSubmesh )
+  {
+    updateTriangulation( _shape );
+
+    Bnd_Box bb;
+    BRepBndLib::Add (_shape, bb);
+    double x1,y1,z1,x2,y2,z2;
+    bb.Get (x1,y1,z1,x2,y2,z2);
+    netgen::Point<3> p1 = netgen::Point<3> (x1,y1,z1);
+    netgen::Point<3> p2 = netgen::Point<3> (x2,y2,z2);
+    occgeo.boundingbox = netgen::Box<3> (p1,p2);
+
+    occgeo.shape    = _shape;
+    occgeo.changed  = 1;
+    if ( !useFMapFunction )
+    {
+      int totNbFaces = 0;
+      TopAbs_ShapeEnum types[4] = { TopAbs_VERTEX, TopAbs_EDGE, TopAbs_FACE, TopAbs_SOLID };
+      for ( TopAbs_ShapeEnum shType : types )
       {
-        if ( !internals.hasInternalVertexInFace() ) {
-          FillSMesh( occgeo, *_ngMesh, initState, *_mesh, nodeVec, comment );
-          initState = NETGENPlugin_ngMeshInfo(_ngMesh);
-        }
-        SMESH_ProxyMesh::Ptr viscousMesh;
-        SMESH_MesherHelper   helper( *_mesh );
-        for ( int faceID = 1; faceID <= occgeo.fmap.Extent(); ++faceID )
+        TopTools_IndexedMapOfShape shapes;
+        TopExp::MapShapes( _shape, shType, shapes );
+        if ( shType == TopAbs_FACE )
+          totNbFaces += 1;
+        for ( int i = 1; i <= shapes.Size(); ++i )
         {
-          const TopoDS_Face& F = TopoDS::Face( occgeo.fmap( faceID ));
-          viscousMesh = StdMeshers_ViscousLayers2D::Compute( *_mesh, F );
-          if ( !viscousMesh )
-            return false;
-          if ( viscousMesh->NbProxySubMeshes() == 0 )
-            continue;
-          // exclude from computation ng segments built on EDGEs of F
-          for (int i = 1; i <= _ngMesh->GetNSeg(); i++)
-          {
-            netgen::Segment & seg = _ngMesh->LineSegment(i);
-            if (seg.si == faceID)
-              seg.si = 0;
+          TopoDS_Shape shape = shapes(i);
+          if ( shape.Orientation() >= TopAbs_INTERNAL )
+            shape.Orientation( TopAbs_FORWARD );
+
+          switch ( shType ) {
+            case TopAbs_FACE  : occgeo.fmap.Add( shape ); break;
+            case TopAbs_EDGE  : occgeo.emap.Add( shape ); break;
+            case TopAbs_VERTEX: occgeo.vmap.Add( shape ); break;
+            case TopAbs_SOLID :occgeo.somap.Add( shape ); break;
+            default:;
           }
-          // add new segments to _ngMesh instead of excluded ones
-          helper.SetSubShape( F );
-          TSideVector wires =
-            StdMeshers_FaceSide::GetFaceWires( F, *_mesh, /*skipMediumNodes=*/true,
-                                               error, &helper, viscousMesh );
-          error = AddSegmentsToMesh( *_ngMesh, occgeo, wires, helper, nodeVec );
-
-          if ( !error ) error = SMESH_ComputeError::New();
         }
-        initState = NETGENPlugin_ngMeshInfo(_ngMesh);
-      }
-
-      // Let netgen compute 2D mesh
-      startWith = netgen::MESHCONST_MESHSURFACE;
-      endWith = _optimize ? netgen::MESHCONST_OPTSURFACE : netgen::MESHCONST_MESHSURFACE;
-      try
-      {
-        OCC_CATCH_SIGNALS;
-
-        err = ngLib.GenerateMesh(occgeo, startWith, endWith);
-
-        if ( netgen::multithread.terminate )
-          return false;
-
-        comment << text (err);
-      }
-      catch (Standard_Failure& ex)
-      {
-        comment << text(ex);
-        //err = 1; -- try to make volumes anyway
-      }
-      catch (netgen::NgException& exc)
-      {
-        comment << text(exc);
-        //err = 1; -- try to make volumes anyway
       }
+      occgeo.facemeshstatus.SetSize (totNbFaces);
+      occgeo.facemeshstatus = 0;
+      occgeo.face_maxh_modified.SetSize(totNbFaces);
+      occgeo.face_maxh_modified = 0;
+      occgeo.face_maxh.SetSize(totNbFaces);
     }
-    if ( _isVolume )
+    else
     {
-      doneTime += faceMeshingTime + ( _optimize ? faceOptimizTime : 0 );
-      _ticTime = doneTime / _totalTime / _progressTic;
+      occgeo.BuildFMap(); // It is possible to let netgen::OCCGeometry to auto fill the maps, it generates slightly different results for the automeshing
     }
-    // ---------------------
-    // generate volume mesh
-    // ---------------------
-    // Fill _ngMesh with nodes and faces of computed 2D submeshes
-    if ( !err && _isVolume &&
-         ( !meshedSM[ MeshDim_2D ].empty() || mparams.quad || _viscousLayersHyp ))
-    {
-      // load SMESH with computed segments and faces
-      FillSMesh( occgeo, *_ngMesh, initState, *_mesh, nodeVec, comment, &quadHelper );
+  }
+  // else
+  // {
+  //   std::cout << "Branch of PrepareOCCgeometry\n";
+  //   PrepareOCCgeometry(occgeo, _mesh->GetShapeToMesh(), *_mesh, meshedSM, internals );
+  // }
+    
+  occgeo.face_maxh = netgen::mparam.maxh;  
+  _occgeom = &occgeo;
 
-      // compute prismatic boundary volumes
-      smIdType nbQuad = _mesh->NbQuadrangles();
-      SMESH_ProxyMesh::Ptr viscousMesh;
-      if ( _viscousLayersHyp )
-      {
-        viscousMesh = _viscousLayersHyp->Compute( *_mesh, _shape );
-        if ( !viscousMesh )
-          return false;
-      }
-      // compute pyramids on quadrangles
-      vector<SMESH_ProxyMesh::Ptr> pyramidMeshes( occgeo.somap.Extent() );
-      if ( nbQuad > 0 )
-        for ( int iS = 1; iS <= occgeo.somap.Extent(); ++iS )
-        {
-          StdMeshers_QuadToTriaAdaptor* adaptor = new StdMeshers_QuadToTriaAdaptor;
-          pyramidMeshes[ iS-1 ].reset( adaptor );
-          bool ok = adaptor->Compute( *_mesh, occgeo.somap(iS), viscousMesh.get() );
-          if ( !ok )
-            return false;
-        }
-      // add proxy faces to NG mesh
-      list< SMESH_subMesh* > viscousSM;
-      for ( int iS = 1; iS <= occgeo.somap.Extent(); ++iS )
-      {
-        list< SMESH_subMesh* > quadFaceSM;
-        for (TopExp_Explorer face(occgeo.somap(iS), TopAbs_FACE); face.More(); face.Next())
-          if ( pyramidMeshes[iS-1] && pyramidMeshes[iS-1]->GetProxySubMesh( face.Current() ))
-          {
-            quadFaceSM.push_back( _mesh->GetSubMesh( face.Current() ));
-            meshedSM[ MeshDim_2D ].remove( quadFaceSM.back() );
-          }
-          else if ( viscousMesh && viscousMesh->GetProxySubMesh( face.Current() ))
-          {
-            viscousSM.push_back( _mesh->GetSubMesh( face.Current() ));
-            meshedSM[ MeshDim_2D ].remove( viscousSM.back() );
-          }
-        if ( !quadFaceSM.empty() )
-          FillNgMesh(occgeo, *_ngMesh, nodeVec, quadFaceSM, &quadHelper, pyramidMeshes[iS-1]);
-      }
-      if ( !viscousSM.empty() )
-        FillNgMesh(occgeo, *_ngMesh, nodeVec, viscousSM, &quadHelper, viscousMesh );
+  // Fill basic mesh param and 
+  // define _ngMesh
+  SetBasicMeshParameters( ngLib, mparams, occgeo );
+  // Mesh internal faces and edges with netgen and then fill smesh with those entries!
+  if ( internals->hasInternalEdges() )
+    FillInternalElements( ngLib, *internals, occgeo, initState, quadHelper, meshedSM );   
+}
 
-      // fill _ngMesh with faces of sub-meshes
-      err = ! ( FillNgMesh(occgeo, *_ngMesh, nodeVec, meshedSM[ MeshDim_2D ], &quadHelper));
-      initState = NETGENPlugin_ngMeshInfo(_ngMesh, /*checkRemovedElems=*/true);
-      // toPython( _ngMesh )
-    }
-    if (!err && _isVolume)
-    {
-      // Pass 3D simple parameters to NETGEN
-      const NETGENPlugin_SimpleHypothesis_3D* simple3d =
-        dynamic_cast< const NETGENPlugin_SimpleHypothesis_3D* > ( _simpleHyp );
-      if ( simple3d ) {
-        _ngMesh->Compress();
-        if ( double vol = simple3d->GetMaxElementVolume() ) {
-          // max volume
-          mparams.maxh = pow( 72, 1/6. ) * pow( vol, 1/3. );
-          mparams.maxh = min( mparams.maxh, occgeo.boundingbox.Diam()/2 );
-        }
-        else {
-          // length from faces
-          mparams.maxh = _ngMesh->AverageH();
-        }
-        _ngMesh->SetGlobalH (mparams.maxh);
-        mparams.grading = 0.4;
-        ngLib.CalcLocalH( ngLib._ngMesh );
-      }
-      // Care of vertices internal in solids and internal faces (issue 0020676)
-      if ( internals.hasInternalVertexInSolid() || internals.hasInternalFaces() )
+void NETGENPlugin_Mesher::FillSMESH( netgen::OCCGeometry& occgeo, NETGENPlugin_ngMeshInfo& initState, 
+                                      vector< const SMDS_MeshNode* >& nodeVec, SMESH_MesherHelper &quadHelper,
+                                      SMESH_Comment& comment )
+{
+  FillSMesh( occgeo, *_ngMesh, initState, *_mesh, nodeVec, comment, &quadHelper );
+
+  if ( quadHelper.GetIsQuadratic() ) // remove free nodes
+  {
+    for ( size_t i = 0; i < nodeVec.size(); ++i )
+      if ( nodeVec[i] && nodeVec[i]->NbInverseElements() == 0 )
       {
-        // store computed faces in SMESH in order not to create SMESH
-        // faces for ng faces added here
-        FillSMesh( occgeo, *_ngMesh, initState, *_mesh, nodeVec, comment, &quadHelper );
-        // add ng faces to solids with internal vertices
-        AddIntVerticesInSolids( occgeo, *_ngMesh, nodeVec, internals );
-        // duplicate mesh faces on internal faces
-        FixIntFaces( occgeo, *_ngMesh, internals );
-        initState = NETGENPlugin_ngMeshInfo(_ngMesh);
+        _mesh->GetMeshDS()->RemoveFreeNode( nodeVec[i], 0, /*fromGroups=*/false );
+        nodeVec[i]=0;
       }
-      // Let netgen compute 3D mesh
-      startWith = endWith = netgen::MESHCONST_MESHVOLUME;
-      try
-      {
-        OCC_CATCH_SIGNALS;
+    for ( size_t i = nodeVec.size()-1; i > 0; --i ) // remove trailing removed nodes
+      if ( !nodeVec[i] )
+        nodeVec.resize( i );
+      else
+        break;
+  }
+}
 
-        err = ngLib.GenerateMesh(occgeo, startWith, endWith);
+bool NETGENPlugin_Mesher::Compute( NETGENPlugin_NetgenLibWrapper& ngLib, vector< const SMDS_MeshNode* >& nodeVec, bool write2SMESH, DIM dim )
+{
+  netgen::OCCGeometry occgeo;
+  list< SMESH_subMesh* > meshedSM[3]; 
+  NETGENPlugin_Internals internals( *_mesh, _shape, _isVolume );  
+  NETGENPlugin_ngMeshInfo initState; // it remembers size of ng mesh equal to size of Smesh
+  netgen::MeshingParameters& mparams = netgen::mparam;  
+  SMESH_MesherHelper quadHelper( *_mesh );
+  quadHelper.SetIsQuadratic( mparams.secondorder );
+  SMESH_ComputeErrorPtr error = SMESH_ComputeError::New();
+  _ngMesh = NULL;
+  SMESH_Comment comment;
 
-        if ( netgen::multithread.terminate )
-          return false;
+  InitialSetupSA( ngLib, occgeo, meshedSM, &internals, quadHelper, initState, mparams, true );
+  int err = Fill0D1DElements( occgeo, nodeVec, meshedSM, quadHelper );  
+  initState = NETGENPlugin_ngMeshInfo(_ngMesh);
 
-        if ( comment.empty() ) // do not overwrite a previous error
-          comment << text(err);
-      }
-      catch (Standard_Failure& ex)
-      {
-        if ( comment.empty() ) // do not overwrite a previous error
-          comment << text(ex);
-        err = 1;
-      }
-      catch (netgen::NgException& exc)
-      {
-        if ( comment.empty() ) // do not overwrite a previous error
-          comment << text(exc);
-        err = 1;
-      }
-      _ticTime = ( doneTime += voluMeshingTime ) / _totalTime / _progressTic;
+  if ( dim >= 1 )
+    err = Compute1D( ngLib, occgeo );
+  if( dim >= 2 )
+    err = Compute2D( ngLib, occgeo, mparams, meshedSM, initState, &internals, nodeVec, comment, dim );
+  if ( dim == 3 )
+    err = Compute3D( ngLib, occgeo, mparams, meshedSM, initState, &internals, nodeVec, quadHelper, comment );
 
-      // Let netgen optimize 3D mesh
-      if ( !err && _optimize )
-      {
-        startWith = endWith = netgen::MESHCONST_OPTVOLUME;
-        try
-        {
-          OCC_CATCH_SIGNALS;
+  if ( write2SMESH )
+  {
+    _mesh->Clear();
+    FillSMESH( occgeo, initState, nodeVec, quadHelper, comment );
+  }
 
-          err = ngLib.GenerateMesh(occgeo, startWith, endWith);
+  return err;
+}
 
-          if ( netgen::multithread.terminate )
-            return false;
+bool NETGENPlugin_Mesher::Compute1D( NETGENPlugin_NetgenLibWrapper& ngLib, netgen::OCCGeometry& occgeo )
+{
+  return CallNetgenMeshEdges( ngLib, occgeo );
+}
 
-          if ( comment.empty() ) // do not overwrite a previous error
-            comment << text(err);
-        }
-        catch (Standard_Failure& ex)
-        {
-          if ( comment.empty() ) // do not overwrite a previous error
-            comment << text(ex);
-        }
-        catch (netgen::NgException& exc)
-        {
-          if ( comment.empty() ) // do not overwrite a previous error
-            comment << text(exc);
-        }
-      }
-    }
-    if (!err && mparams.secondorder > 0)
+bool NETGENPlugin_Mesher::Compute2D( NETGENPlugin_NetgenLibWrapper& ngLib, netgen::OCCGeometry& occgeo,
+                                      netgen::MeshingParameters &mparams, list< SMESH_subMesh* >* meshedSM, 
+                                      NETGENPlugin_ngMeshInfo& initState, NETGENPlugin_Internals* internals, 
+                                      vector< const SMDS_MeshNode* >& nodeVec, SMESH_Comment& comment, DIM dim )
+{                                           
+  SetBasicMeshParametersFor2D(occgeo, nodeVec, mparams, internals, initState );
+  if ( !Fill2DViscousLayer(occgeo, nodeVec, internals, initState ) )
+    return false;
+
+  mparams.uselocalh = true; // restore as it is used at surface optimization
+  int err = CallNetgenMeshFaces( ngLib, occgeo, comment );
+  
+  if ( !err && dim == DIM::D2 /* if mesh is 3D then second order is defined after volumens are computed*/ )
+    MakeSecondOrder( ngLib, mparams, occgeo, meshedSM, initState, comment );
+
+  return err;
+}
+
+bool NETGENPlugin_Mesher::Compute3D( NETGENPlugin_NetgenLibWrapper& ngLib, netgen::OCCGeometry& occgeo,
+                                    netgen::MeshingParameters &mparams, list< SMESH_subMesh* >* meshedSM, 
+                                    NETGENPlugin_ngMeshInfo& initState, NETGENPlugin_Internals* internals, 
+                                    vector< const SMDS_MeshNode* >& nodeVec, SMESH_MesherHelper &quadHelper, 
+                                    SMESH_Comment& comment )
+{
+  int err = 0;
+  if ( !Fill3DViscousLayerAndQuadAdaptor( occgeo, nodeVec, mparams, initState, meshedSM, quadHelper, err ) )
+    return false;
+
+  if ( !err )
+  {
+    SetBasicMeshParametersFor3D( ngLib, occgeo, nodeVec, mparams, internals, initState, quadHelper, comment );
+    err = CallNetgenMeshVolumens( ngLib, occgeo, comment );
+  }
+
+  if ( !err )
+    MakeSecondOrder( ngLib, mparams, occgeo, meshedSM, initState, comment );
+
+  return err;
+}
+
+//=============================================================================
+/*!
+ * Here we are going to use the NETGEN mesher
+ */
+//=============================================================================
+
+bool NETGENPlugin_Mesher::Compute()
+{
+  NETGENPlugin_NetgenLibWrapper ngLib;
+  netgen::OCCGeometry occgeo;
+  list< SMESH_subMesh* > meshedSM[3]; // for 0-2 dimensions
+  NETGENPlugin_Internals internals( *_mesh, _shape, _isVolume );  
+  NETGENPlugin_ngMeshInfo initState; // it remembers size of ng mesh equal to size of Smesh
+
+  netgen::MeshingParameters& mparams = netgen::mparam;  
+  SMESH_MesherHelper quadHelper( *_mesh );
+  quadHelper.SetIsQuadratic( mparams.secondorder );
+
+  vector< const SMDS_MeshNode* > nodeVec;
+  SMESH_ComputeErrorPtr error = SMESH_ComputeError::New();
+  _ngMesh = NULL;
+  
+  int err = 0;
+
+  _totalTime = edgeFaceMeshingTime;
+  if ( _optimize )
+    _totalTime += faceOptimizTime;
+  if ( _isVolume )
+    _totalTime += voluMeshingTime + ( _optimize ? volOptimizeTime : 0 );
+  double doneTime = 0;
+  _ticTime = -1;
+  _progressTic = 1;
+  _curShapeIndex = -1;
+
+  // -------------------------
+  // Generate the mesh
+  // -------------------------
+
+  InitialSetup( ngLib, occgeo, meshedSM, &internals, quadHelper, initState, mparams );
+  err = Fill0D1DElements( occgeo, nodeVec, meshedSM, quadHelper );  
+  initState = NETGENPlugin_ngMeshInfo(_ngMesh);
+  err = CallNetgenMeshEdges( ngLib, occgeo );
+  SetBasicMeshParametersFor2D(occgeo, nodeVec, mparams, &internals, initState );
+  if ( !Fill2DViscousLayer(occgeo, nodeVec, &internals, initState ) )
+    return false;
+
+  SMESH_Comment comment;
+  
+  {
+
+    if ( _isVolume )
+      _ticTime = ( doneTime += edgeMeshingTime ) / _totalTime / _progressTic;
+
+    mparams.uselocalh = true; // restore as it is used at surface optimization
+    err = CallNetgenMeshFaces( ngLib, occgeo, comment );
+    
+    if ( _isVolume )
     {
-      try
-      {
-        OCC_CATCH_SIGNALS;
-        if ( !meshedSM[ MeshDim_1D ].empty() )
-        {
-          // remove segments not attached to geometry (IPAL0052479)
-          for (int i = 1; i <= _ngMesh->GetNSeg(); ++i)
-          {
-            const netgen::Segment & seg = _ngMesh->LineSegment (i);
-            if ( seg.epgeominfo[ 0 ].edgenr == 0 )
-            {
-              _ngMesh->DeleteSegment( i );
-              initState._nbSegments--;
-            }
-          }
-          _ngMesh->Compress();
-        }
-        // convert to quadratic
-#ifdef NETGEN_V6
-        occgeo.GetRefinement().MakeSecondOrder(*_ngMesh);
-#else
-        netgen::OCCRefinementSurfaces(occgeo).MakeSecondOrder(*_ngMesh);
-#endif
+      doneTime += faceMeshingTime + ( _optimize ? faceOptimizTime : 0 );
+      _ticTime = doneTime / _totalTime / _progressTic;
+    }
 
-        // care of elements already loaded to SMESH
-        // if ( initState._nbSegments > 0 )
-        //   makeQuadratic( occgeo.emap, _mesh );
-        // if ( initState._nbFaces > 0 )
-        //   makeQuadratic( occgeo.fmap, _mesh );
-      }
-      catch (Standard_Failure& ex)
-      {
-        if ( comment.empty() ) // do not overwrite a previous error
-          comment << "Exception in netgen at passing to 2nd order ";
-      }
-      catch (netgen::NgException& exc)
-      {
-        if ( comment.empty() ) // do not overwrite a previous error
-          comment << exc.What();
-      }
+    if ( !err )
+      if ( !Fill3DViscousLayerAndQuadAdaptor( occgeo, nodeVec, mparams, initState, meshedSM, quadHelper, err ) )
+        return false;
+
+    if (!err && _isVolume)
+    {
+      SetBasicMeshParametersFor3D( ngLib, occgeo, nodeVec, mparams, &internals, initState, quadHelper, comment );
+      err = CallNetgenMeshVolumens( ngLib, occgeo, comment );
     }
+
+    if (!err )
+      MakeSecondOrder( ngLib, mparams, occgeo, meshedSM, initState, comment );
   }
 
   _ticTime = 0.98 / _progressTic;
@@ -3414,23 +3671,9 @@ bool NETGENPlugin_Mesher::Compute()
   // Feed back the SMESHDS with the generated Nodes and Elements
   if ( true /*isOK*/ ) // get whatever built
   {
-    FillSMesh( occgeo, *_ngMesh, initState, *_mesh, nodeVec, comment, &quadHelper );
-
-    if ( quadHelper.GetIsQuadratic() ) // remove free nodes
-    {
-      for ( size_t i = 0; i < nodeVec.size(); ++i )
-        if ( nodeVec[i] && nodeVec[i]->NbInverseElements() == 0 )
-        {
-          _mesh->GetMeshDS()->RemoveFreeNode( nodeVec[i], 0, /*fromGroups=*/false );
-          nodeVec[i]=0;
-        }
-      for ( size_t i = nodeVec.size()-1; i > 0; --i ) // remove trailing removed nodes
-        if ( !nodeVec[i] )
-          nodeVec.resize( i );
-        else
-          break;
-    }
+    FillSMESH( occgeo, initState, nodeVec, quadHelper, comment );    
   }
+
   SMESH_ComputeErrorPtr readErr = ReadErrors(nodeVec);
   if ( readErr && readErr->HasBadElems() )
   {
@@ -4097,7 +4340,6 @@ NETGENPlugin_Internals::NETGENPlugin_Internals( SMESH_Mesh&         mesh,
   for ( f.Init( shape, TopAbs_FACE ); f.More(); f.Next() )
   {
     int faceID = meshDS->ShapeToIndex( f.Current() );
-
     // find not computed internal edges
 
     for ( e.Init( f.Current().Oriented(TopAbs_FORWARD), TopAbs_EDGE ); e.More(); e.Next() )
index fe66ad2623ed5a3a17f1ee3e4617996995ad776c..f33b154f63d52eff1d4fd031ac92b57a43124a52 100644 (file)
 #include "SALOME_Basics.hxx"
 #include <TopTools_IndexedMapOfShape.hxx>
 
+// Netgen include files
+#ifndef OCCGEOMETRY
+#define OCCGEOMETRY
+#endif
+#include <occgeom.hpp>
+#include <meshing.hpp>
+
 namespace nglib {
 #include <nglib.h>
 }
@@ -60,6 +67,8 @@ class TopoDS_Shape;
 namespace netgen {
   class OCCGeometry;
   class Mesh;
+  NETGENPLUGIN_DLL_HEADER
+  extern MeshingParameters mparam;
 }
 
 // Class for temporary folder switching
@@ -172,7 +181,56 @@ class NETGENPLUGIN_EXPORT NETGENPlugin_Mesher
   void SetLocalSizeForChordalError( netgen::OCCGeometry& occgeo, netgen::Mesh& ngMesh );
   static void SetLocalSize( netgen::OCCGeometry& occgeo, netgen::Mesh& ngMesh );
 
+  
+/**
+ * @brief InitialSetup. Fill occgeo map with geometrical objects not meshed. Fill meshdSM with the already computed
+ *        submeshes, and mesh the internal edges so faces with internal are eventurally properly meshed.
+ *        Define the class members _ngMesh and _occgeom
+ */
+  void InitialSetup( NETGENPlugin_NetgenLibWrapper& ngLib, netgen::OCCGeometry& occgeo, 
+                      list< SMESH_subMesh* >* meshedSM, NETGENPlugin_Internals* internals, 
+                      SMESH_MesherHelper &quadHelper, NETGENPlugin_ngMeshInfo& initState, netgen::MeshingParameters &mparams );
+
+  void InitialSetupSA( NETGENPlugin_NetgenLibWrapper& ngLib, netgen::OCCGeometry& occgeo, 
+                      list< SMESH_subMesh* >* meshedSM, NETGENPlugin_Internals* internals, 
+                      SMESH_MesherHelper &quadHelper, NETGENPlugin_ngMeshInfo& initState, 
+                      netgen::MeshingParameters &mparams, bool useFMapFunction = false );
+
+  void SetBasicMeshParameters( NETGENPlugin_NetgenLibWrapper& ngLib, netgen::MeshingParameters &mparams, netgen::OCCGeometry& occgeo );
+  void SetBasicMeshParametersFor2D( netgen::OCCGeometry& occgeo, vector< const SMDS_MeshNode* >& nodeVec, 
+                                      netgen::MeshingParameters &mparams, NETGENPlugin_Internals* internals, 
+                                        NETGENPlugin_ngMeshInfo& initState );
+  void SetBasicMeshParametersFor3D( NETGENPlugin_NetgenLibWrapper& ngLib, netgen::OCCGeometry& occgeo, 
+                                      vector< const SMDS_MeshNode* >& nodeVec, netgen::MeshingParameters &mparams, 
+                                        NETGENPlugin_Internals* internals, NETGENPlugin_ngMeshInfo& initState, SMESH_MesherHelper &quadHelper, 
+                                          SMESH_Comment& comment );
+  void CallNetgenConstAnalysis( NETGENPlugin_NetgenLibWrapper& ngLib, netgen::MeshingParameters &mparams, netgen::OCCGeometry& occgeo );
+  int CallNetgenMeshEdges( NETGENPlugin_NetgenLibWrapper& ngLib, netgen::OCCGeometry& occgeo );
+  int CallNetgenMeshFaces( NETGENPlugin_NetgenLibWrapper& ngLib, netgen::OCCGeometry& occgeo, SMESH_Comment& comment );
+  int CallNetgenMeshVolumens( NETGENPlugin_NetgenLibWrapper& ngLib, netgen::OCCGeometry& occgeo, SMESH_Comment& comment );
+  void MakeSecondOrder( NETGENPlugin_NetgenLibWrapper& ngLib, netgen::MeshingParameters &mparams, netgen::OCCGeometry& occgeo, 
+                          list< SMESH_subMesh* >* meshedSM, NETGENPlugin_ngMeshInfo& initState, SMESH_Comment& comment );
+  int FillInternalElements( NETGENPlugin_NetgenLibWrapper& ngLib, NETGENPlugin_Internals& internals, netgen::OCCGeometry& occgeo,
+                              NETGENPlugin_ngMeshInfo& initState, SMESH_MesherHelper &quadHelper, list< SMESH_subMesh* >* meshedSM );
+  bool Fill2DViscousLayer( netgen::OCCGeometry& occgeo, vector< const SMDS_MeshNode* >& nodeVec, 
+                            NETGENPlugin_Internals* internals, NETGENPlugin_ngMeshInfo& initState );
+
+  bool Fill3DViscousLayerAndQuadAdaptor( netgen::OCCGeometry& occgeo, vector< const SMDS_MeshNode* >& nodeVec, 
+                                          netgen::MeshingParameters &mparams, NETGENPlugin_ngMeshInfo& initState,
+                                          list< SMESH_subMesh* >* meshedSM, SMESH_MesherHelper &quadHelper, int & err );
+
+  int Fill0D1DElements( netgen::OCCGeometry& occgeo, vector< const SMDS_MeshNode* >& nodeVec, list< SMESH_subMesh* >* meshedSM, SMESH_MesherHelper &quadHelper );
+  void FillSMESH( netgen::OCCGeometry& occgeo, NETGENPlugin_ngMeshInfo& initState, vector< const SMDS_MeshNode* >& nodeVec, SMESH_MesherHelper &quadHelper, SMESH_Comment& comment );
+  ///// End definition methods to rewrite function
+  
+  enum DIM {
+    D1 = 1,
+    D2,
+    D3
+  };       
+  
   bool Compute();
+  bool Compute( NETGENPlugin_NetgenLibWrapper& ngLib, vector< const SMDS_MeshNode* >& nodeVec, bool write2SMESH, DIM dim );
 
   bool Evaluate(MapShapeNbElems& aResMap);
 
@@ -244,6 +302,19 @@ class NETGENPLUGIN_EXPORT NETGENPlugin_Mesher
 
  private:
 
+  bool Compute1D( NETGENPlugin_NetgenLibWrapper& ngLib, netgen::OCCGeometry& occgeo );
+
+  bool Compute2D( NETGENPlugin_NetgenLibWrapper& ngLib, netgen::OCCGeometry& occgeo,
+                  netgen::MeshingParameters &mparams, list< SMESH_subMesh* >* meshedSM, 
+                  NETGENPlugin_ngMeshInfo& initState, NETGENPlugin_Internals* internals, 
+                  vector< const SMDS_MeshNode* >& nodeVec, SMESH_Comment& comment, DIM dim );
+                  
+  bool Compute3D( NETGENPlugin_NetgenLibWrapper& ngLib, netgen::OCCGeometry& occgeo,
+                  netgen::MeshingParameters &mparams, list< SMESH_subMesh* >* meshedSM, 
+                  NETGENPlugin_ngMeshInfo& initState, NETGENPlugin_Internals* internals, 
+                  vector< const SMDS_MeshNode* >& nodeVec, SMESH_MesherHelper &quadHelper, 
+                  SMESH_Comment& comment);
+
   SMESH_Mesh*          _mesh;
   const TopoDS_Shape&  _shape;
   bool                 _isVolume;
diff --git a/src/NETGENPlugin/NETGENPlugin_NETGEN_1D2D3D_SA.cxx b/src/NETGENPlugin/NETGENPlugin_NETGEN_1D2D3D_SA.cxx
new file mode 100644 (file)
index 0000000..bcaa38c
--- /dev/null
@@ -0,0 +1,357 @@
+// Copyright (C) 2007-2023  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
+//
+
+// NETGENPlugin : C++ implementation
+// File   : NETGENPlugin_NETGEN_2D3D_SA.cxx
+// Author : Cesar Conopoima (OCC)
+// Date   : 01/11/2023
+// Project   : SALOME
+//=============================================================================
+//
+
+#include "NETGENPlugin_DriverParam.hxx"
+#include "NETGENPlugin_Hypothesis.hxx"
+#include "NETGENPlugin_SimpleHypothesis_3D.hxx"
+
+#include "NETGENPlugin_NETGEN_2D.hxx"
+#include "NETGENPlugin_NETGEN_1D2D3D_SA.hxx"
+
+#include <SMESHDS_Mesh.hxx>
+#include <SMESH_ControlsDef.hxx>
+#include <SMESH_Gen.hxx>
+#include <SMESH_Mesh.hxx>
+#include <SMESH_MesherHelper.hxx>
+#include <SMESH_DriverShape.hxx>
+#include <SMESH_DriverMesh.hxx>
+#include <utilities.h>
+
+#include <list>
+
+namespace nglib {
+#include <nglib.h>
+}
+#include <meshing.hpp>
+
+using namespace nglib;
+
+//=============================================================================
+/*!
+ *  
+ */
+//=============================================================================
+
+NETGENPlugin_NETGEN_1D2D3D_SA::NETGENPlugin_NETGEN_1D2D3D_SA()
+  : NETGENPlugin_NETGEN_2D3D(0, new SMESH_Gen())
+{
+  _name = "NETGEN_1D2D3D_SA";  
+}
+
+//=============================================================================
+/*!
+ *  
+ */
+//=============================================================================
+
+NETGENPlugin_NETGEN_1D2D3D_SA::~NETGENPlugin_NETGEN_1D2D3D_SA()
+{
+}
+
+/**
+ * @brief fill plugin hypothesis from the netgen_params structure
+ *
+ * @param aParams the structure
+ */
+void NETGENPlugin_NETGEN_1D2D3D_SA::fillHyp(netgen_params aParams)
+{  
+ if( aParams.myType == hypoType::Hypo )
+  {
+    NETGENPlugin_Hypothesis * hypParameters = new NETGENPlugin_Hypothesis(0, GetGen());
+
+    hypParameters->SetMaxSize(aParams.maxh);
+    hypParameters->SetMinSize(aParams.minh);
+    hypParameters->SetNbSegPerEdge(aParams.segmentsperedge);
+    hypParameters->SetGrowthRate(aParams.grading);
+    hypParameters->SetNbSegPerRadius(aParams.curvaturesafety);
+    hypParameters->SetSecondOrder(aParams.secondorder);
+    hypParameters->SetQuadAllowed(aParams.quad);
+    hypParameters->SetOptimize(aParams.optimize);
+    hypParameters->SetFineness((NETGENPlugin_Hypothesis::Fineness)aParams.fineness);
+    hypParameters->SetSurfaceCurvature(aParams.uselocalh);
+    hypParameters->SetFuseEdges(aParams.merge_solids);
+    hypParameters->SetChordalErrorEnabled(aParams.chordalError);
+    if(aParams.optimize){
+      hypParameters->SetNbSurfOptSteps(aParams.optsteps2d);
+      hypParameters->SetNbVolOptSteps(aParams.optsteps3d);
+    }
+    hypParameters->SetElemSizeWeight(aParams.elsizeweight);
+    hypParameters->SetWorstElemMeasure(aParams.opterrpow);
+    hypParameters->SetUseDelauney(aParams.delaunay);
+    hypParameters->SetCheckOverlapping(aParams.checkoverlap);
+    hypParameters->SetCheckChartBoundary(aParams.checkchartboundary);
+    hypParameters->SetMeshSizeFile(aParams.meshsizefilename);
+
+    _hypothesis = dynamic_cast< const NETGENPlugin_Hypothesis *> (hypParameters);
+  }
+  else if ( aParams.myType == hypoType::Simple2D )
+  {
+    NETGENPlugin_SimpleHypothesis_2D * hypParameters = new NETGENPlugin_SimpleHypothesis_2D(0, GetGen());
+    
+    // mandatory to fill in this branch case!
+    // Number of segments   (int)
+    // localLenght          (double)
+    // maxElement area      (double)
+    // GetAllowQuadrangles  (bool)
+    hypParameters->SetNumberOfSegments( aParams.numberOfSegments );
+    if ( !aParams.numberOfSegments )
+      hypParameters->SetLocalLength( aParams.localLength );      
+    hypParameters->SetMaxElementArea( aParams.maxElementArea );
+    hypParameters->SetAllowQuadrangles( aParams.allowQuadrangles );
+    _hypothesis = dynamic_cast< const NETGENPlugin_SimpleHypothesis_2D *> (hypParameters);
+  }
+  else if ( aParams.myType == hypoType::Simple3D )
+  {
+    NETGENPlugin_SimpleHypothesis_3D * hypParameters = new NETGENPlugin_SimpleHypothesis_3D(0, GetGen());
+    
+    // mandatory to fill in this branch case!
+    // Number of segments   (int)
+    // localLenght          (double)
+    // maxElement area      (double)
+    // maxElement volume    (double)
+    // GetAllowQuadrangles  (bool)
+    hypParameters->SetNumberOfSegments( aParams.numberOfSegments );
+    if ( !aParams.numberOfSegments )
+      hypParameters->SetLocalLength( aParams.localLength );      
+    hypParameters->SetMaxElementArea( aParams.maxElementArea );
+    hypParameters->SetMaxElementVolume( aParams.maxElementVol );
+    hypParameters->SetAllowQuadrangles( aParams.allowQuadrangles );
+
+    _hypothesis = dynamic_cast< const NETGENPlugin_SimpleHypothesis_3D *> (hypParameters);
+  }
+}
+
+/**
+ * @brief Write a binary file containing information on the elements/nodes
+ *        created by the netgen mesher
+ *
+ * @param nodeVec mapping between the mesh id and the netgen structure id
+ * @param ngLib Wrapper on netgen library
+ * @param new_element_file Name of the output file
+ * @param NumOfPremeshedNodes Number of nodes in the netgen structure
+ * @return true if there are some error
+ */
+bool NETGENPlugin_NETGEN_1D2D3D_SA::FillNewElementFile( std::vector< const SMDS_MeshNode* > &nodeVec, NETGENPlugin_NetgenLibWrapper &ngLib, 
+                                                        std::string new_element_file, const NETGENPlugin_Mesher::DIM dim )
+{
+
+  // Particularities: As submeshing is not supported nodeVect is empty and NumberOfPremeshedNodes is also zero
+  Ng_Mesh* NetgenMesh  = ngLib.ngMesh();
+  int NetgenNodes       = Ng_GetNP(NetgenMesh);
+  int NetgenSeg2D       = Ng_GetNSeg_2D( NetgenMesh );
+  int NetgenFaces       = Ng_GetNSE(NetgenMesh);
+  int NetgenVols        = Ng_GetNE(NetgenMesh);
+  int segmentId;
+  bool isOK = ( NetgenNodes > 0 && NetgenSeg2D > 0 );
+  if ( isOK && !new_element_file.empty() )
+  {
+    MESSAGE("Writting new elements")
+    std::ofstream df(new_element_file, ios::out|ios::binary);
+
+    double NetgenPoint[3];
+    int    NetgenSegment[2];
+    int    NetgenSurface[8];
+    int    NetgenVolumens[10];
+
+    // Writing nodevec (correspondence netgen numbering mesh numbering)
+    // Number of nodes
+    const int NumOfPremeshedNodes = nodeVec.size();
+    df.write((char*) &NumOfPremeshedNodes, sizeof(int));
+    df.write((char*) &NetgenNodes, sizeof(int));
+
+    for (int nodeIndex = 1 ; nodeIndex <= NumOfPremeshedNodes; ++nodeIndex )
+    {
+      //Id of the point
+      int id = nodeVec.at(nodeIndex)->GetID();
+      df.write((char*) &id, sizeof(int));
+    }
+
+    // Writing all new points
+    for (int nodeIndex = NumOfPremeshedNodes + 1; nodeIndex <= NetgenNodes; ++nodeIndex )
+    {
+      Ng_GetPoint( NetgenMesh, nodeIndex, NetgenPoint );
+      // Coordinates of the point
+      df.write((char *) &NetgenPoint, sizeof(double)*3);
+    }
+
+    if ( dim >= NETGENPlugin_Mesher::D1 )
+    {
+      // create segments at boundaries.
+      df.write((char*) &NetgenSeg2D, sizeof(int));
+      for ( int elemIndex = 1; elemIndex <= NetgenSeg2D; ++elemIndex )
+      {
+        Ng_GetSegment_2D( NetgenMesh, elemIndex, NetgenSegment, &segmentId );
+        df.write((char*) &NetgenSegment, sizeof(int) * 2 );
+      }
+    }
+    if ( dim >= NETGENPlugin_Mesher::D2 ) 
+    {
+      // create surface elements.
+      df.write((char*) &NetgenFaces, sizeof(int));
+      for ( int elemIndex = 1; elemIndex <= NetgenFaces; ++elemIndex )
+      {
+        nglib::Ng_Surface_Element_Type elemType = Ng_GetSurfaceElement( NetgenMesh, elemIndex, NetgenSurface );
+        switch (elemType)
+        {
+          case nglib::NG_TRIG:
+          {
+            df.write((char*) &NetgenSurface, sizeof(int) * 3 );
+            break;
+          }            
+          case nglib::NG_QUAD:
+          {
+            df.write((char*) &NetgenSurface, sizeof(int) * 4 );
+            break;
+          }            
+          case nglib::NG_TRIG6:
+          {
+            df.write((char*) &NetgenSurface, sizeof(int) * 6 );
+            break;
+          }            
+          case nglib::NG_QUAD8:
+          {
+            df.write((char*) &NetgenSurface, sizeof(int) * 8 );
+            break;
+          }
+          default:
+          { break; }         
+        }
+      }
+    }
+    if ( dim >= NETGENPlugin_Mesher::D3 ) 
+    {
+      // create volume elements.
+      df.write((char*) &NetgenVols, sizeof(int));
+      for ( int elemIndex = 1; elemIndex <= NetgenVols; ++elemIndex )
+      {
+        nglib::Ng_Volume_Element_Type elemType =  Ng_GetVolumeElement( NetgenMesh, elemIndex, NetgenVolumens );
+        switch (elemType)
+        {        
+          case nglib::NG_TET:
+          {
+            df.write((char*) &NetgenVolumens, sizeof(int) * 4 );
+            break;
+          } 
+          case nglib::NG_PYRAMID:
+          {
+            df.write((char*) &NetgenVolumens, sizeof(int) * 5 );
+            break;
+          }
+          case nglib::NG_PRISM:
+          {
+            df.write((char*) &NetgenVolumens, sizeof(int) * 6 );
+            break;
+          }
+          case nglib::NG_TET10:
+          {
+            df.write((char*) &NetgenVolumens, sizeof(int) * 10 );
+            break;
+          }
+          default:
+          { break; }
+        }
+      }
+    }          
+    df.close();
+  }
+  return false;
+}
+
+/**
+ * @brief Compute the mesh based on the 
+ *
+ * @param aMesh the read Mesh
+ * @param aShape the loaded shape
+ * @param new_element_file output file containing info the elements created by the mesher
+ * @param output_mesh whether or not write the created elements into the mesh
+ * @param dim the dimension to be meshed.
+ * @return negation of mesh fail: true, false
+ */
+bool NETGENPlugin_NETGEN_1D2D3D_SA::Compute(SMESH_Mesh& aMesh, TopoDS_Shape &aShape, std::string new_element_file, bool output_mesh, NETGENPlugin_Mesher::DIM dim )
+{
+  //
+  netgen::multithread.terminate = 0;
+  NETGENPlugin_Mesher mesher(&aMesh, aShape, /*is3D = */ false );
+  mesher.SetParameters(dynamic_cast<const NETGENPlugin_Hypothesis*>(_hypothesis));
+  if ( dim == NETGENPlugin_Mesher::D3 )
+    mesher.SetParameters(dynamic_cast<const NETGENPlugin_SimpleHypothesis_3D*>(_hypothesis));
+  else
+    mesher.SetParameters(dynamic_cast<const NETGENPlugin_SimpleHypothesis_2D*>(_hypothesis));
+  NETGENPlugin_NetgenLibWrapper ngLib;
+  vector< const SMDS_MeshNode* > nodeVec;
+  bool err = mesher.Compute( ngLib, nodeVec, output_mesh, dim );  
+  FillNewElementFile( nodeVec, ngLib, new_element_file, dim );
+  return err;
+}
+
+/**
+ * @brief Running the mesher on the given files
+ *
+ * @param input_mesh_file Mesh file (containing 2D elements)
+ * @param shape_file Shape file (BREP or STEP format)
+ * @param hypo_file Ascii file containing the netgen parameters
+ * @param element_orientation_file Binary file containing the orientation of surface elemnts
+ * @param new_element_file output file containing info the elements created by the mesher
+ * @param output_mesh_file output mesh file (if empty it will not be created)
+ * @return negation of mesh fail: true, false
+ */
+int NETGENPlugin_NETGEN_1D2D3D_SA::run(const std::string input_mesh_file,
+                                      const std::string shape_file,
+                                      const std::string hypo_file,
+                                      const std::string element_orientation_file,
+                                      const std::string new_element_file,
+                                      const std::string output_mesh_file, 
+                                      const NETGENPlugin_Mesher::DIM dim )
+  {
+
+  // _element_orientation_file = element_orientation_file;
+
+  std::unique_ptr<SMESH_Mesh> myMesh(_gen->CreateMesh(false));
+
+  // Importing mesh
+  SMESH_DriverMesh::importMesh(input_mesh_file, *myMesh);
+  // Importing shape
+  TopoDS_Shape myShape;
+  SMESH_DriverShape::importShape(shape_file, myShape);
+  // Importing hypothesis
+  netgen_params myParams;  
+  importNetgenParams(hypo_file, myParams);
+  fillHyp(myParams);
+
+  int ret = Compute( *myMesh, myShape, new_element_file, !output_mesh_file.empty(), dim );
+  if(ret){
+    std::cerr << "Meshing failed" << std::endl;
+    return ret;
+  }
+
+  if(!output_mesh_file.empty()){
+    std::string meshName = "MESH";
+    SMESH_DriverMesh::exportMesh(output_mesh_file, *myMesh, meshName);
+  }
+
+  return ret;
+}
diff --git a/src/NETGENPlugin/NETGENPlugin_NETGEN_1D2D3D_SA.hxx b/src/NETGENPlugin/NETGENPlugin_NETGEN_1D2D3D_SA.hxx
new file mode 100644 (file)
index 0000000..c269567
--- /dev/null
@@ -0,0 +1,62 @@
+// Copyright (C) 2007-2023  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
+//
+
+// NETGENPlugin : C++ implementation
+// File      : NETGENPlugin_NETGEN_2D_SA.hxx
+// Author    : Cesar Conopoima (OCC)
+// Date      : 23/10/2023
+// Project   : SALOME
+//=============================================================================
+//
+#ifndef _NETGENPlugin_NETGEN_1D2D3D_SA_HXX_
+#define _NETGENPlugin_NETGEN_1D2D3D_SA_HXX_
+
+#include "NETGENPlugin_Defs.hxx"
+#include "NETGENPlugin_Mesher.hxx"
+#include "NETGENPlugin_NETGEN_2D3D.hxx"
+#include "SMESH_Algo.hxx"
+#include "SMESH_Mesh.hxx"
+
+class NETGENPlugin_Mesher;
+
+class NETGENPLUGIN_EXPORT NETGENPlugin_NETGEN_1D2D3D_SA: public NETGENPlugin_NETGEN_2D3D
+{
+public:
+  NETGENPlugin_NETGEN_1D2D3D_SA();
+  virtual ~NETGENPlugin_NETGEN_1D2D3D_SA();
+
+  bool Compute(SMESH_Mesh& aMesh, TopoDS_Shape &aShape, std::string new_element_file, bool output_mesh, NETGENPlugin_Mesher::DIM dim );
+  
+  int run(const std::string input_mesh_file,
+          const std::string shape_file,
+          const std::string hypo_file,
+          const std::string element_orientation_file,
+          const std::string new_element_file,
+          const std::string output_mesh_file,
+          const NETGENPlugin_Mesher::DIM dim );
+
+protected:
+  
+  void fillHyp(netgen_params aParams);
+  bool FillNewElementFile( std::vector< const SMDS_MeshNode* > &nodeVec, NETGENPlugin_NetgenLibWrapper &ngLib,
+                            std::string new_element_file, const NETGENPlugin_Mesher::DIM dim );
+
+};
+
+#endif
index 16213e86aa9a684bfe210312eae33dc5754e5d2f..4c4e1e6ab20cf51c17c301cee5163f2de2aa04ce 100644 (file)
@@ -22,8 +22,6 @@
 // Project   : SALOME
 //
 #include "NETGENPlugin_NETGEN_2D_ONLY.hxx"
-
-#include "NETGENPlugin_Mesher.hxx"
 #include "NETGENPlugin_Hypothesis_2D.hxx"
 
 #include <SMDS_MeshElement.hxx>
@@ -46,6 +44,8 @@
 
 #include <utilities.h>
 
+#include <GEOMUtils.hxx>
+
 #include <list>
 #include <vector>
 #include <limits>
@@ -61,6 +61,23 @@ namespace nglib {
 #endif
 #include <occgeom.hpp>
 #include <meshing.hpp>
+
+#include "SMESH_Octree.hxx"
+
+//// Used for node projection in curve
+#include <BRepAdaptor_Curve.hxx>
+#include <BRep_Builder.hxx>
+#include <BRep_Tool.hxx>
+#include <BndLib_Add3dCurve.hxx>
+#include <GCPnts_TangentialDeflection.hxx>
+#include <ShapeAnalysis_Curve.hxx>
+#include <TopExp.hxx>
+#include <TopExp_Explorer.hxx>
+#include <TopoDS.hxx>
+#include <TopoDS_Compound.hxx>
+#include <TopoDS_Edge.hxx>
+#include <TopoDS_Vertex.hxx>
+
 //#include <meshtype.hpp>
 namespace netgen {
   NETGENPLUGIN_DLL_HEADER
@@ -74,6 +91,264 @@ using namespace std;
 using namespace netgen;
 using namespace nglib;
 
+namespace // copied class from StdMesher_Importe_1D   
+{
+  /*!
+  * \brief Compute point position on a curve. Use octree to fast reject far points
+  */
+  class CurveProjector : public SMESH_Octree
+  {
+  public:
+    CurveProjector( const TopoDS_Edge& edge, double enlarge );
+
+    bool IsOnCurve( const gp_XYZ& point, double & distance2, double & u );
+
+    bool IsOut( const gp_XYZ& point ) const { return getBox()->IsOut( point ); }
+
+  protected:
+    CurveProjector() {}
+    SMESH_Octree* newChild() const { return new CurveProjector; }
+    void          buildChildrenData();
+    Bnd_B3d*      buildRootBox();
+
+  private:
+    struct CurveSegment : public Bnd_B3d
+    {
+      double _chord, _chord2, _length2;
+      gp_Pnt _pFirst, _pLast;
+      gp_Lin _line;
+      Handle(Geom_Curve) _curve;
+
+      CurveSegment() {}
+      void Init( const gp_Pnt& pf, const gp_Pnt& pl,
+                 double uf, double ul, double tol, Handle(Geom_Curve)& curve );
+      bool IsOn( const gp_XYZ& point, double & distance2, double & u );
+      bool IsInContact( const Bnd_B3d& bb );
+    };
+    std::vector< CurveSegment > _segments;
+  };
+
+    //===============================================================================
+  /*!
+   * \brief Create an octree of curve segments
+   */
+  //================================================================================
+
+  CurveProjector::CurveProjector( const TopoDS_Edge& edge, double enlarge )
+    :SMESH_Octree( 0 )
+  {
+    double f,l;
+    Handle(Geom_Curve) curve = BRep_Tool::Curve( edge, f, l );
+    double curDeflect = 0.3; // Curvature deflection
+    double angDeflect = 1e+100; // Angular deflection - don't control chordal error
+    GCPnts_TangentialDeflection div( BRepAdaptor_Curve( edge ), angDeflect, curDeflect );
+    _segments.resize( div.NbPoints() - 1 );
+    for ( int i = 1; i < div.NbPoints(); ++i )
+      try {
+        _segments[ i - 1 ].Init( div.Value( i ),     div.Value( i+1 ),
+                                 div.Parameter( i ), div.Parameter( i+1 ),
+                                 enlarge, curve );
+      }
+      catch ( Standard_Failure ) {
+        _segments.resize( _segments.size() - 1 );
+        --i;
+      }
+    if ( _segments.size() < 3 )
+      myIsLeaf = true;
+
+    compute();
+
+    if ( _segments.size() == 1 )
+      myBox->Enlarge( enlarge );
+  }
+
+  //================================================================================
+  /*!
+   * \brief Return the maximal box
+   */
+  //================================================================================
+
+  Bnd_B3d* CurveProjector::buildRootBox()
+  {
+    Bnd_B3d* box = new Bnd_B3d;
+    for ( size_t i = 0; i < _segments.size(); ++i )
+      box->Add( _segments[i] );
+    return box;
+  }
+
+  //================================================================================
+  /*!
+   * \brief Redistribute segments among children
+   */
+  //================================================================================
+
+  void CurveProjector::buildChildrenData()
+  {
+    bool allIn = true;
+    for ( size_t i = 0; i < _segments.size(); ++i )
+    {
+      for (int j = 0; j < 8; j++)
+      {
+        if ( _segments[i].IsInContact( *myChildren[j]->getBox() ))
+          ((CurveProjector*)myChildren[j])->_segments.push_back( _segments[i]);
+        else
+          allIn = false;
+      }
+    }
+    if ( allIn && _segments.size() < 3 )
+    {
+      myIsLeaf = true;
+      for (int j = 0; j < 8; j++)
+        static_cast<CurveProjector*>( myChildren[j])->myIsLeaf = true;
+    }
+    else
+    {
+      SMESHUtils::FreeVector( _segments ); // = _segments.clear() + free memory
+
+      for (int j = 0; j < 8; j++)
+      {
+        CurveProjector* child = static_cast<CurveProjector*>( myChildren[j]);
+        if ( child->_segments.size() < 3 )
+          child->myIsLeaf = true;
+      }
+    }
+  }
+
+  //================================================================================
+  /*!
+   * \brief Return true if a point is close to the curve
+   *  \param [in] point - the point
+   *  \param [out] distance2 - distance to the curve
+   *  \param [out] u - parameter on the curve
+   *  \return bool - is the point is close to the curve
+   */
+  //================================================================================
+
+  bool CurveProjector::IsOnCurve( const gp_XYZ& point, double & distance2, double & u )
+  {
+    if ( getBox()->IsOut( point ))
+      return false;
+
+    bool ok = false;
+    double dist2, param;
+    distance2 = Precision::Infinite();
+
+    if ( isLeaf() )
+    {
+      for ( size_t i = 0; i < _segments.size(); ++i )
+        if ( !_segments[i].IsOut( point ) &&
+             _segments[i].IsOn( point, dist2, param ) &&
+             dist2 < distance2 )
+        {
+          distance2 = dist2;
+          u         = param;
+          ok        = true;
+        }
+      return ok;
+    }
+    else
+    {
+      for (int i = 0; i < 8; i++)
+        if (((CurveProjector*) myChildren[i])->IsOnCurve( point, dist2, param ) &&
+            dist2 < distance2 )
+        {
+          distance2 = dist2;
+          u         = param;
+          ok        = true;
+        }
+    }
+    return ok;
+  }
+
+  //================================================================================
+  /*!
+   * \brief Initialize
+   */
+  //================================================================================
+
+  void CurveProjector::CurveSegment::Init(const gp_Pnt&       pf,
+                                          const gp_Pnt&       pl,
+                                          const double        uf,
+                                          const double        ul,
+                                          const double        tol,
+                                          Handle(Geom_Curve)& curve )
+  {
+    _pFirst  = pf;
+    _pLast   = pl;
+    _curve   = curve;
+    _length2 = pf.SquareDistance( pl );
+    _line.SetLocation( pf );
+    _line.SetDirection( gp_Vec( pf, pl ));
+    _chord2  = Max( _line.     SquareDistance( curve->Value( uf + 0.25 * ( ul - uf ))),
+                    Max( _line.SquareDistance( curve->Value( uf + 0.5  * ( ul - uf ))),
+                         _line.SquareDistance( curve->Value( uf + 0.75 * ( ul - uf )))));
+    _chord2 *= ( 1.05 * 1.05 ); // +5%
+    _chord2  = Max( tol, _chord2 );
+    _chord   = Sqrt( _chord2 );
+
+    Bnd_Box bb;
+    BndLib_Add3dCurve::Add( GeomAdaptor_Curve( curve, uf, ul ), tol, bb );
+    Add( bb.CornerMin() );
+    Add( bb.CornerMax() );
+  }
+
+  //================================================================================
+  /*!
+   * \brief Return true if a point is close to the curve segment
+   *  \param [in] point - the point
+   *  \param [out] distance2 - distance to the curve
+   *  \param [out] u - parameter on the curve
+   *  \return bool - is the point is close to the curve segment
+   */
+  //================================================================================
+
+  bool CurveProjector::CurveSegment::IsOn( const gp_XYZ& point, double & distance2, double & u )
+  {
+    distance2 = _line.SquareDistance( point );
+    if ( distance2 > _chord2 )
+      return false;
+
+    // check if the point projection falls into the segment range
+    {
+      gp_Vec edge( _pFirst, _pLast );
+      gp_Vec n1p ( _pFirst, point  );
+      u = ( edge * n1p ) / _length2; // param [0,1] on the edge
+      if ( u < 0. )
+      {
+        if ( _pFirst.SquareDistance( point ) > _chord2 )
+          return false;
+      }
+      else if ( u > 1. )
+      {
+        if ( _pLast.SquareDistance( point ) > _chord2 )
+          return false;
+      }
+    }
+    gp_Pnt proj;
+    distance2 = ShapeAnalysis_Curve().Project( _curve, point, Precision::Confusion(),
+                                               proj, u, false );
+    distance2 *= distance2;
+    return true;
+  }
+
+  //================================================================================
+  /*!
+   * \brief Check if the segment is in contact with a box
+   */
+  //================================================================================
+
+  bool CurveProjector::CurveSegment::IsInContact( const Bnd_B3d& bb )
+  {
+    if ( bb.IsOut( _line.Position(), /*isRay=*/true, _chord ))
+      return false;
+
+    gp_Ax1 axRev = _line.Position().Reversed();
+    axRev.SetLocation( _pLast );
+    return !bb.IsOut( axRev, /*isRay=*/true, _chord );
+  }
+
+}
+
 //=============================================================================
 /*!
  *
@@ -217,43 +492,216 @@ bool NETGENPlugin_NETGEN_2D_ONLY::CheckHypothesis (SMESH_Mesh&         aMesh,
 //   }
 // }
 
-//=============================================================================
-/*!
- *Here we are going to use the NETGEN mesher
- */
-//=============================================================================
 
-bool NETGENPlugin_NETGEN_2D_ONLY::Compute(SMESH_Mesh&         aMesh,
-                                          const TopoDS_Shape& aShape)
-{
-  netgen::multithread.terminate = 0;
-  //netgen::multithread.task = "Surface meshing";
 
-  SMESHDS_Mesh* meshDS = aMesh.GetMeshDS();
+/**
+ * @brief MapSegmentsToEdges. 
+ * @remark To feed 1D segments not associated to any geometry we need:
+ *          1) For each face, check all segments that are in the face (use ShapeAnalysis_Surface class)
+ *          2) Check to which edge the segment below to, use the copied [from StdMesher_Importe_1D] CurveProjector class 
+ *          3) Create new netgen segments with the (u,v) parameters obtained from the ShapeAnalysis_Surface projector
+ *          4) also define the 'param' value of the nodes relative to the edges obtained from CurveProjector
+ *          5) Add the new netgen segments IN ORDER into the netgen::mesh data structure to form a closed chain
+ *          6) Beware with the occ::edge orientation
+ * 
+ * @param aMesh Mesh file (containing 1D elements)
+ * @param aShape Shape file (BREP or STEP format)
+ * @param ngLib netgenlib library wrapper
+ * @param nodeVec vector of nodes used internally to feed smesh aMesh after computation
+ * @param premeshedNodes map of prmeshed nodes and the smesh nodeID associate to it
+ * @param newNetgenCoordinates map of 3D coordinate of new points created by netgen
+ * @param newNetgenElements map of triangular or quadrangular elements ID and the nodes defining the 2D element
+ * @return false
+ */
+bool NETGENPlugin_NETGEN_2D_ONLY::MapSegmentsToEdges(SMESH_Mesh& aMesh, const TopoDS_Shape& aShape, NETGENPlugin_NetgenLibWrapper &ngLib,
+                                                      vector< const SMDS_MeshNode* >& nodeVec, std::map<int,const SMDS_MeshNode*>& premeshedNodes, 
+                                                      std::map<int,std::vector<double>>& newNetgenCoordinates, std::map<int,std::vector<smIdType>>& newNetgenElements )
+{  
+  SMESHDS_Mesh * meshDS = aMesh.GetMeshDS();
   SMESH_MesherHelper helper(aMesh);
   helper.SetElementsOnShape( true );
+  const int numberOfPremeshedNodes = aMesh.NbNodes();
+  TopTools_IndexedMapOfShape faces;
+  TopExp::MapShapes( aShape, TopAbs_FACE, faces );
 
-  NETGENPlugin_NetgenLibWrapper ngLib;
-  ngLib._isComputeOk = false;
+  for ( int i = 1; i <= faces.Size(); ++i )
+  {
+    int numOfEdges = 0;
+    int totalEdgeLenght = 0;
 
-  netgen::Mesh   ngMeshNoLocSize;
-  netgen::Mesh * ngMeshes[2] = { (netgen::Mesh*) ngLib._ngMesh,  & ngMeshNoLocSize };
-  netgen::OCCGeometry occgeoComm;
+    TopoDS_Face face = TopoDS::Face(faces( i ) );
+    Bnd_Box FaceBox;
+    GEOMUtils::PreciseBoundingBox( face, FaceBox );
+    if ( face.Orientation() != TopAbs_FORWARD && face.Orientation() != TopAbs_REVERSED )
+      face.Orientation( TopAbs_FORWARD );
 
-  // min / max sizes are set as follows:
-  // if ( _hypParameters )
-  //    min and max are defined by the user
-  // else if ( _hypLengthFromEdges )
-  //    min = aMesher.GetDefaultMinSize()
-  //    max = average segment len of a FACE
-  // else if ( _hypMaxElementArea )
-  //    min = aMesher.GetDefaultMinSize()
-  //    max = f( _hypMaxElementArea )
-  // else
-  //    min = aMesher.GetDefaultMinSize()
-  //    max = max segment len of a FACE
+    // Set the occgeom to be meshed and some ngMesh parameteres
+    netgen::OCCGeometry occgeom;
+    netgen::Mesh * ngMesh =  (netgen::Mesh*) ngLib._ngMesh;
+    ngMesh->DeleteMesh();
 
-  NETGENPlugin_Mesher aMesher( &aMesh, aShape, /*isVolume=*/false);
+    occgeom.shape = face;
+    occgeom.fmap.Add( face );
+    occgeom.CalcBoundingBox();
+    occgeom.facemeshstatus.SetSize(1);
+    occgeom.facemeshstatus = 0;
+    occgeom.face_maxh_modified.SetSize(1);
+    occgeom.face_maxh_modified = 0;
+    occgeom.face_maxh.SetSize(1);
+    occgeom.face_maxh = netgen::mparam.maxh;
+
+    // Set the face descriptor
+    const int solidID = 0, faceID = 1; /*always 1 because faces are meshed one by one*/
+    if ( ngMesh->GetNFD() < 1 )
+      ngMesh->AddFaceDescriptor( netgen::FaceDescriptor( faceID, solidID, solidID, 0 ));
+
+    ngMesh->SetGlobalH ( netgen::mparam.maxh );
+    ngMesh->SetMinimalH( netgen::mparam.minh );
+    Box<3> bb = occgeom.GetBoundingBox();
+    bb.Increase (bb.Diam()/10);
+    ngMesh->SetLocalH (bb.PMin(), bb.PMax(), netgen::mparam.grading);  
+    // end set the occgeom to be meshed and some ngMesh parameteres
+
+    Handle(ShapeAnalysis_Surface) sprojector = new ShapeAnalysis_Surface( BRep_Tool::Surface( face ));
+    double tol = BRep_Tool::MaxTolerance( face, TopAbs_FACE );
+    gp_Pnt surfPnt(0,0,0);
+  
+    map<const SMDS_MeshNode*, int > node2ngID;
+    map<int, const SMDS_MeshNode* > ng2smesh;
+    // Figure out which edge is this onde in!
+    TopTools_IndexedMapOfShape edges;
+    TopExp::MapShapes( face, TopAbs_EDGE, edges );
+    TopoDS_Edge meshingEdge;
+    // Check wich nodes are in this face!
+    SMDS_ElemIteratorPtr iteratorElem = meshDS->elementsIterator(SMDSAbs_Edge);
+    std::unique_ptr<CurveProjector> myCurveProjector;
+    while ( iteratorElem->more() ) // loop on elements on a geom face
+    {
+      // check mesh face
+      const SMDS_MeshElement* elem = iteratorElem->next();
+      const SMDS_MeshNode* node0 = elem->GetNode( 0 );
+      const SMDS_MeshNode* node1 = elem->GetNode( 1 );
+      SMESH_NodeXYZ nXYZ0( node0 );
+      SMESH_NodeXYZ nXYZ1( node1 );
+      double segmentLength = ( nXYZ0 - nXYZ1 ).Modulus();
+      nXYZ0 += nXYZ1;
+      nXYZ0 /= 2.0;
+      if( FaceBox.IsOut( nXYZ0 ) )
+        continue;      
+
+      gp_XY uv = sprojector->ValueOfUV( nXYZ0, tol ).XY();
+      surfPnt = sprojector->Value( uv );
+      double dist = surfPnt.Distance( nXYZ0 );
+      
+      // used for the curve projector of edges
+      double geomTol = Precision::Confusion();
+      if ( dist < tol /*element is on face*/ )
+      {
+        numOfEdges++;
+        totalEdgeLenght += segmentLength;
+        int occEdgeIdFound = -1;
+        for ( int edgeId = 1; edgeId <= edges.Size(); ++edgeId ) /*find in which edge the node is placed*/
+        {
+          meshingEdge = TopoDS::Edge(edges( edgeId ));
+          myCurveProjector = std::unique_ptr<CurveProjector>( new CurveProjector(meshingEdge, geomTol) );
+          if ( myCurveProjector->IsOut( nXYZ0 ) /*keep searching*/)
+            continue;
+          else
+          {
+            occEdgeIdFound = edgeId;
+            break;
+          }
+        }
+
+        netgen::Segment seg;
+        for ( size_t nodeId = 0; nodeId < 2; nodeId++)
+        {
+          const SMDS_MeshNode* node = elem->GetNode( nodeId );
+          int ngId = ngMesh->GetNP() + 1;
+          ngId = node2ngID.insert( make_pair( node, ngId )).first->second;
+          if ( ngId > ngMesh->GetNP() /* mean it is a new node to be add to the mesh*/)
+          {
+            // Restric size of mesh based on the edge dimension
+            {
+              SMESH_NodeXYZ nXYZ( node );
+              netgen::Point3d pi(nXYZ.X(), nXYZ.Y(), nXYZ.Z());
+              ngMesh->RestrictLocalH( pi, segmentLength );
+            }
+              
+            netgen::MeshPoint mp( netgen::Point<3> (node->X(), node->Y(), node->Z()) );
+            ngMesh->AddPoint ( mp, 1, netgen::EDGEPOINT );
+            ng2smesh.insert( make_pair(ngId, node) );
+            premeshedNodes.insert( make_pair( (int)node->GetID(), node ) );
+          } 
+          seg[nodeId] = ngId;                // ng node id
+          SMESH_NodeXYZ nXYZ( node );
+          gp_XY uv    = sprojector->ValueOfUV( nXYZ, tol ).XY();
+
+          // Compute the param (relative distance) of the node in relation the edge
+          // fundamental for netgen working properly
+          double dist2, param;
+          if ( myCurveProjector->IsOnCurve( nXYZ, dist2, param ) )
+          {
+            seg.epgeominfo[ nodeId ].dist = param;
+            seg.epgeominfo[ nodeId ].u    = uv.X();
+            seg.epgeominfo[ nodeId ].v    = uv.Y();   
+            seg.epgeominfo[ nodeId ].edgenr = occEdgeIdFound;
+          }       
+        }
+
+        if ( meshingEdge.Orientation() != TopAbs_FORWARD )
+        {
+          swap(seg[0], seg[1]);
+          swap(seg.epgeominfo[0], seg.epgeominfo[1] );
+        }
+        
+        seg.edgenr = ngMesh->GetNSeg() + 1; // netgen segment id
+        seg.si     = faceID;
+        ngMesh->AddSegment( seg );
+      } // end if edge is on the face    
+    } // end iteration on elements
+    
+    // set parameters from _hypLengthFromEdges if needed
+    if ( !_hypParameters && _hypLengthFromEdges && numOfEdges )
+    {
+      netgen::mparam.maxh = totalEdgeLenght / numOfEdges;
+      if ( netgen::mparam.maxh < DBL_MIN )
+        netgen::mparam.maxh = occgeom.GetBoundingBox().Diam();
+
+      occgeom.face_maxh = netgen::mparam.maxh;
+      ngMesh->SetGlobalH ( netgen::mparam.maxh );
+    }    
+    // end set parameters
+
+    ngMesh->CalcSurfacesOfNode();
+    const int startWith = MESHCONST_MESHSURFACE;
+    const int endWith   = MESHCONST_OPTSURFACE;
+    int err = ngLib.GenerateMesh(occgeom, startWith, endWith, ngMesh);
+
+    // Ng_Mesh * ngMeshptr = (Ng_Mesh*) ngLib._ngMesh;
+    // int NetgenNodes = Ng_GetNP(ngMeshptr);
+    // int NetgenSeg_2D = Ng_GetNSeg_2D(ngMeshptr);
+    // int NetgenFaces = Ng_GetNSE(ngMeshptr);
+    // std::cout << "\n";
+    // std::cout << "Number of nodes, seg, faces, vols: " << NetgenNodes << ", " << NetgenSeg_2D << ", " << NetgenFaces << "\n";
+    // std::cout << "err from netgen computation: " << err << "\n";
+
+    // ----------------------------------------------------
+    // Fill the SMESHDS with the generated nodes and faces
+    // ----------------------------------------------------
+    nodeVec.clear();
+    FillNodesAndElements( aMesh, helper, ngMesh, nodeVec, ng2smesh, newNetgenCoordinates, newNetgenElements, numberOfPremeshedNodes );
+  } // Face iteration
+
+  return false;
+}
+
+std::tuple<bool,bool> NETGENPlugin_NETGEN_2D_ONLY::SetParameteres( SMESH_Mesh& aMesh, const TopoDS_Shape& aShape, 
+                                                                  NETGENPlugin_Mesher& aMesher, netgen::Mesh * ngMeshes,
+                                                                  netgen::OCCGeometry& occgeoComm, bool isSubMeshSupported )
+{
+  SMESHDS_Mesh* meshDS = aMesh.GetMeshDS();
+    
   aMesher.SetParameters( _hypParameters ); // _hypParameters -> netgen::mparam
   const bool toOptimize = _hypParameters ? _hypParameters->GetOptimize() : true;
   if ( _hypMaxElementArea )
@@ -269,7 +717,6 @@ bool NETGENPlugin_NETGEN_2D_ONLY::Compute(SMESH_Mesh&         aMesh,
 
   if ( isCommonLocalSize ) // compute common local size in ngMeshes[0]
   {
-    //list< SMESH_subMesh* > meshedSM[4]; --> all sub-shapes are added to occgeoComm
     aMesher.PrepareOCCgeometry( occgeoComm, aShape, aMesh );//, meshedSM );
 
     // local size set at MESHCONST_ANALYSE step depends on
@@ -291,44 +738,268 @@ bool NETGENPlugin_NETGEN_2D_ONLY::Compute(SMESH_Mesh&         aMesh,
     occgeoComm.face_maxh = netgen::mparam.maxh;
 #ifdef NETGEN_V6
     netgen::OCCParameters occparam;
-    netgen::OCCSetLocalMeshSize( occgeoComm, *ngMeshes[0], netgen::mparam, occparam );
+    netgen::OCCSetLocalMeshSize( occgeoComm, *ngMeshes, netgen::mparam, occparam );
 #else
-    netgen::OCCSetLocalMeshSize( occgeoComm, *ngMeshes[0] );
+    netgen::OCCSetLocalMeshSize( occgeoComm, *ngMeshes );
 #endif
     occgeoComm.emap.Clear();
     occgeoComm.vmap.Clear();
 
-    // set local size according to size of existing segments
-    TopTools_IndexedMapOfShape edgeMap;
-    TopExp::MapShapes( aMesh.GetShapeToMesh(), TopAbs_EDGE, edgeMap );
-    for ( int iE = 1; iE <= edgeMap.Extent(); ++iE )
+    if ( isSubMeshSupported )
     {
-      const TopoDS_Shape& edge = edgeMap( iE );
-      if ( SMESH_Algo::isDegenerated( TopoDS::Edge( edge )))
-        continue;
-      SMESHDS_SubMesh* smDS = meshDS->MeshElements( edge );
-      if ( !smDS ) continue;
-      SMDS_ElemIteratorPtr segIt = smDS->GetElements();
-      while ( segIt->more() )
+      // set local size according to size of existing segments
+      TopTools_IndexedMapOfShape edgeMap;
+      TopExp::MapShapes( aMesh.GetShapeToMesh(), TopAbs_EDGE, edgeMap );
+      for ( int iE = 1; iE <= edgeMap.Extent(); ++iE )
       {
-        const SMDS_MeshElement* seg = segIt->next();
-        SMESH_TNodeXYZ n1 = seg->GetNode(0);
-        SMESH_TNodeXYZ n2 = seg->GetNode(1);
-        gp_XYZ p = 0.5 * ( n1 + n2 );
+        const TopoDS_Shape& edge = edgeMap( iE );
+        if ( SMESH_Algo::isDegenerated( TopoDS::Edge( edge )))
+          continue;
+        SMESHDS_SubMesh* smDS = meshDS->MeshElements( edge );
+        if ( !smDS ) continue;
+        SMDS_ElemIteratorPtr segIt = smDS->GetElements();
+        while ( segIt->more() )
+        {
+          const SMDS_MeshElement* seg = segIt->next();
+          SMESH_TNodeXYZ n1 = seg->GetNode(0);
+          SMESH_TNodeXYZ n2 = seg->GetNode(1);
+          gp_XYZ p = 0.5 * ( n1 + n2 );
+          netgen::Point3d pi(p.X(), p.Y(), p.Z());
+          ngMeshes->RestrictLocalH( pi, factor * ( n1 - n2 ).Modulus() );
+        }
+      }
+    }
+    else
+    {
+      SMDS_ElemIteratorPtr iteratorElem = meshDS->elementsIterator(SMDSAbs_Edge);
+      while ( iteratorElem->more() ) // loop on elements on a geom face
+      {
+        const SMDS_MeshElement* elem = iteratorElem->next();      
+        const SMDS_MeshNode* node0 = elem->GetNode( 0 );
+        const SMDS_MeshNode* node1 = elem->GetNode( 1 );
+        SMESH_NodeXYZ nXYZ0( node0 );
+        SMESH_NodeXYZ nXYZ1( node1 );
+        double segmentLength = ( nXYZ0 - nXYZ1 ).Modulus();
+        gp_XYZ p = 0.5 * ( nXYZ0 + nXYZ1 );
         netgen::Point3d pi(p.X(), p.Y(), p.Z());
-        ngMeshes[0]->RestrictLocalH( pi, factor * ( n1 - n2 ).Modulus() );
+        ngMeshes->RestrictLocalH( pi, factor * segmentLength );
       }
     }
 
     // set local size defined on shapes
-    aMesher.SetLocalSize( occgeoComm, *ngMeshes[0] );
-    aMesher.SetLocalSizeForChordalError( occgeoComm, *ngMeshes[0] );
+    aMesher.SetLocalSize( occgeoComm, *ngMeshes );
+    aMesher.SetLocalSizeForChordalError( occgeoComm, *ngMeshes );
     try {
-      ngMeshes[0]->LoadLocalMeshSize( mparam.meshsizefilename );
+      ngMeshes->LoadLocalMeshSize( mparam.meshsizefilename );
     } catch (NgException & ex) {
-      return error( COMPERR_BAD_PARMETERS, ex.What() );
+      throw error( COMPERR_BAD_PARMETERS, ex.What() );
+    }
+  }
+
+  return std::make_tuple( isCommonLocalSize, isDefaultHyp );
+}
+
+bool NETGENPlugin_NETGEN_2D_ONLY::ComputeMaxhOfFace( TopoDS_Face& Face, NETGENPlugin_Mesher& aMesher, TSideVector& wires, 
+                                                      netgen::OCCGeometry& occgeoComm, bool isDefaultHyp, bool isCommonLocalSize )
+{
+  size_t nbWires = wires.size();
+  if ( !_hypParameters )
+  {
+    double edgeLength = 0;
+    if (_hypLengthFromEdges )
+    {
+      // compute edgeLength as an average segment length
+      smIdType nbSegments = 0;
+      for ( size_t iW = 0; iW < nbWires; ++iW )
+      {
+        edgeLength += wires[ iW ]->Length();
+        nbSegments += wires[ iW ]->NbSegments();
+      }
+      if ( nbSegments )
+        edgeLength /= double( nbSegments );
+      netgen::mparam.maxh = edgeLength;
+    }
+    else if ( isDefaultHyp )
+    {
+      // set edgeLength by a longest segment
+      double maxSeg2 = 0;
+      for ( size_t iW = 0; iW < nbWires; ++iW )
+      {
+        const UVPtStructVec& points = wires[ iW ]->GetUVPtStruct();
+        if ( points.empty() )
+          return error( COMPERR_BAD_INPUT_MESH );
+        gp_Pnt pPrev = SMESH_TNodeXYZ( points[0].node );
+        for ( size_t i = 1; i < points.size(); ++i )
+        {
+          gp_Pnt p = SMESH_TNodeXYZ( points[i].node );
+          maxSeg2 = Max( maxSeg2, p.SquareDistance( pPrev ));
+          pPrev = p;
+        }
+      }
+      edgeLength = sqrt( maxSeg2 ) * 1.05;
+      netgen::mparam.maxh = edgeLength;
+    }
+    if ( netgen::mparam.maxh < DBL_MIN )
+      netgen::mparam.maxh = occgeoComm.GetBoundingBox().Diam();
+
+    if ( !isCommonLocalSize )
+    {
+      netgen::mparam.minh = aMesher.GetDefaultMinSize( Face, netgen::mparam.maxh );
+    }
+  }
+  return true;
+}
+
+void NETGENPlugin_NETGEN_2D_ONLY::FillNodesAndElements( SMESH_Mesh& aMesh, SMESH_MesherHelper& helper, netgen::Mesh * ngMesh, 
+                                                        vector< const SMDS_MeshNode* >& nodeVec, map<int, const SMDS_MeshNode* >& ng2smesh, 
+                                                        std::map<int,std::vector<double>>& newNetgenCoordinates, 
+                                                        std::map<int,std::vector<smIdType>>& newNetgenElements, const int numberOfPremeshedNodes )
+{
+  SMESHDS_Mesh* meshDS = aMesh.GetMeshDS();
+
+  int nbNodes = ngMesh->GetNP();
+  int nbFaces = ngMesh->GetNSE();
+  nodeVec.resize( nbNodes + 1, 0 );
+  // map to index local node numeration to global numeration used by Remote mesher to write the final global resulting mesh
+  std::map<smIdType,smIdType> local2globalMap;
+  
+  smIdType myNewCoordinateCounter = newNetgenCoordinates.size() > 0 ? newNetgenCoordinates.rbegin()->first + 1: numberOfPremeshedNodes+1;
+  int myNewFaceCounter = newNetgenElements.size() > 0 ? newNetgenElements.rbegin()->first + 1 : 1;
+  // add nodes
+  for ( int ngID = 1; ngID <= nbNodes; ++ngID )
+  {
+    const MeshPoint& ngPoint = ngMesh->Point( ngID );
+    // Check if ngPoint is not already present because was in the premeshed mesh boundary
+    if ( ng2smesh.count( ngID ) == 0 )
+    {
+      std::vector<double> netgenCoordinates = {ngPoint(0), ngPoint(1), ngPoint(2)};
+      SMDS_MeshNode * node = meshDS->AddNode(ngPoint(0), ngPoint(1), ngPoint(2));
+      nodeVec[ ngID ] = node;
+      newNetgenCoordinates.insert( make_pair( myNewCoordinateCounter, std::move(netgenCoordinates)) );
+      local2globalMap.insert( std::make_pair( node->GetID(), myNewCoordinateCounter ) );
+      myNewCoordinateCounter++;
+    }
+    else
+    {
+      nodeVec[ ngID ] = ng2smesh[ ngID ];
+      local2globalMap.insert( std::make_pair( nodeVec[ ngID ]->GetID(), nodeVec[ ngID ]->GetID() ) );
+    }
+  }
+  // create faces
+  int i,j;
+  vector<const SMDS_MeshNode*> nodes;
+  for ( i = 1; i <= nbFaces ; ++i )
+  {
+    const Element2d& elem = ngMesh->SurfaceElement(i);
+    nodes.resize( elem.GetNP() );
+    for (j=1; j <= elem.GetNP(); ++j)
+    {
+      int pind = elem.PNum(j);
+      if ( pind < 1 )
+        break;
+      nodes[ j-1 ] = nodeVec[ pind ];
+    }
+    if ( j > elem.GetNP() )
+    {
+      std::vector<smIdType> netgenCoordinates = { local2globalMap[nodes[0]->GetID()], local2globalMap[nodes[1]->GetID()], local2globalMap[nodes[2]->GetID()] };
+      newNetgenElements.insert( std::make_pair( myNewFaceCounter, std::move( netgenCoordinates ) ) );
+      helper.AddFace(nodes[0],nodes[1],nodes[2]); 
+      myNewFaceCounter++;    
     }
   }
+}
+
+void NETGENPlugin_NETGEN_2D_ONLY::FillNodesAndElements( SMESH_Mesh& aMesh, SMESH_MesherHelper& helper, netgen::Mesh * ngMesh, vector< const SMDS_MeshNode* >& nodeVec, int faceId )
+{
+  SMESHDS_Mesh* meshDS = aMesh.GetMeshDS();
+
+  int nbNodes = ngMesh->GetNP();
+  int nbFaces = ngMesh->GetNSE();
+
+  int nbInputNodes = (int) nodeVec.size()-1;
+  nodeVec.resize( nbNodes+1, 0 );
+
+  // add nodes
+  for ( int ngID = nbInputNodes + 1; ngID <= nbNodes; ++ngID )
+  {
+    const MeshPoint& ngPoint = ngMesh->Point( ngID );
+    SMDS_MeshNode * node = meshDS->AddNode(ngPoint(0), ngPoint(1), ngPoint(2));
+    nodeVec[ ngID ] = node;
+  }
+
+  // create faces
+  int i,j;
+  vector<const SMDS_MeshNode*> nodes;
+  for ( i = 1; i <= nbFaces ; ++i )
+  {
+    const Element2d& elem = ngMesh->SurfaceElement(i);
+    nodes.resize( elem.GetNP() );
+    for (j=1; j <= elem.GetNP(); ++j)
+    {
+      int pind = elem.PNum(j);
+      if ( pind < 1 )
+        break;
+      nodes[ j-1 ] = nodeVec[ pind ];
+      if ( nodes[ j-1 ]->GetPosition()->GetTypeOfPosition() == SMDS_TOP_3DSPACE )
+      {
+        const PointGeomInfo& pgi = elem.GeomInfoPi(j);
+        meshDS->SetNodeOnFace( nodes[ j-1 ], faceId, pgi.u, pgi.v);
+      }
+    }
+    if ( j > elem.GetNP() )
+    {
+      if ( elem.GetType() == TRIG )
+        helper.AddFace(nodes[0],nodes[1],nodes[2]);
+      else
+        helper.AddFace(nodes[0],nodes[1],nodes[2],nodes[3]);
+    }
+  }
+}
+  
+
+//=============================================================================
+/*!
+ *Here we are going to use the NETGEN mesher
+ */
+//=============================================================================
+
+bool NETGENPlugin_NETGEN_2D_ONLY::Compute(SMESH_Mesh&         aMesh,
+                                          const TopoDS_Shape& aShape)
+{
+  netgen::multithread.terminate = 0;
+  //netgen::multithread.task = "Surface meshing";
+
+  SMESHDS_Mesh* meshDS = aMesh.GetMeshDS();
+  SMESH_MesherHelper helper(aMesh);
+  helper.SetElementsOnShape( true );
+
+  NETGENPlugin_NetgenLibWrapper ngLib;
+  ngLib._isComputeOk = false;
+
+  netgen::Mesh   ngMeshNoLocSize;
+  netgen::Mesh * ngMeshes[2] = { (netgen::Mesh*) ngLib._ngMesh,  & ngMeshNoLocSize };
+
+  netgen::OCCGeometry occgeoComm;
+
+  // min / max sizes are set as follows:
+  // if ( _hypParameters )
+  //    min and max are defined by the user
+  // else if ( _hypLengthFromEdges )
+  //    min = aMesher.GetDefaultMinSize()
+  //    max = average segment len of a FACE
+  // else if ( _hypMaxElementArea )
+  //    min = aMesher.GetDefaultMinSize()
+  //    max = f( _hypMaxElementArea )
+  // else
+  //    min = aMesher.GetDefaultMinSize()
+  //    max = max segment len of a FACE
+
+  NETGENPlugin_Mesher aMesher( &aMesh, aShape, /*isVolume=*/false);
+  auto options = SetParameteres( aMesh, aShape, aMesher, ngMeshes[0], occgeoComm );
+  const bool isCommonLocalSize  = std::get<0>( options );
+  const bool isDefaultHyp       = std::get<1>( options );
+  const bool toOptimize = _hypParameters ? _hypParameters->GetOptimize() : true;
+
   netgen::mparam.uselocalh = toOptimize; // restore as it is used at surface optimization
 
   // ==================
@@ -382,51 +1053,10 @@ bool NETGENPlugin_NETGEN_2D_ONLY::Compute(SMESH_Mesh&         aMesh,
     // compute maxh of a FACE
     // ----------------------
 
-    if ( !_hypParameters )
-    {
-      double edgeLength = 0;
-      if (_hypLengthFromEdges )
-      {
-        // compute edgeLength as an average segment length
-        smIdType nbSegments = 0;
-        for ( size_t iW = 0; iW < nbWires; ++iW )
-        {
-          edgeLength += wires[ iW ]->Length();
-          nbSegments += wires[ iW ]->NbSegments();
-        }
-        if ( nbSegments )
-          edgeLength /= double( nbSegments );
-        netgen::mparam.maxh = edgeLength;
-      }
-      else if ( isDefaultHyp )
-      {
-        // set edgeLength by a longest segment
-        double maxSeg2 = 0;
-        for ( size_t iW = 0; iW < nbWires; ++iW )
-        {
-          const UVPtStructVec& points = wires[ iW ]->GetUVPtStruct();
-          if ( points.empty() )
-            return error( COMPERR_BAD_INPUT_MESH );
-          gp_Pnt pPrev = SMESH_TNodeXYZ( points[0].node );
-          for ( size_t i = 1; i < points.size(); ++i )
-          {
-            gp_Pnt p = SMESH_TNodeXYZ( points[i].node );
-            maxSeg2 = Max( maxSeg2, p.SquareDistance( pPrev ));
-            pPrev = p;
-          }
-        }
-        edgeLength = sqrt( maxSeg2 ) * 1.05;
-        netgen::mparam.maxh = edgeLength;
-      }
-      if ( netgen::mparam.maxh < DBL_MIN )
-        netgen::mparam.maxh = occgeoComm.GetBoundingBox().Diam();
-
-      if ( !isCommonLocalSize )
-      {
-        netgen::mparam.minh = aMesher.GetDefaultMinSize( F, netgen::mparam.maxh );
-      }
-    }
-
+    bool setMaxh = ComputeMaxhOfFace( F, aMesher, wires, occgeoComm, isDefaultHyp, isCommonLocalSize );
+    if (!setMaxh)
+      return setMaxh;
+      
     // prepare occgeom
     netgen::OCCGeometry occgeom;
     occgeom.shape = F;
@@ -546,52 +1176,10 @@ bool NETGENPlugin_NETGEN_2D_ONLY::Compute(SMESH_Mesh&         aMesh,
         }
       }
 
-
       // ----------------------------------------------------
       // Fill the SMESHDS with the generated nodes and faces
       // ----------------------------------------------------
-
-      int nbNodes = ngMesh->GetNP();
-      int nbFaces = ngMesh->GetNSE();
-
-      int nbInputNodes = (int) nodeVec.size()-1;
-      nodeVec.resize( nbNodes+1, 0 );
-
-      // add nodes
-      for ( int ngID = nbInputNodes + 1; ngID <= nbNodes; ++ngID )
-      {
-        const MeshPoint& ngPoint = ngMesh->Point( ngID );
-        SMDS_MeshNode * node = meshDS->AddNode(ngPoint(0), ngPoint(1), ngPoint(2));
-        nodeVec[ ngID ] = node;
-      }
-
-      // create faces
-      int i,j;
-      vector<const SMDS_MeshNode*> nodes;
-      for ( i = 1; i <= nbFaces ; ++i )
-      {
-        const Element2d& elem = ngMesh->SurfaceElement(i);
-        nodes.resize( elem.GetNP() );
-        for (j=1; j <= elem.GetNP(); ++j)
-        {
-          int pind = elem.PNum(j);
-          if ( pind < 1 )
-            break;
-          nodes[ j-1 ] = nodeVec[ pind ];
-          if ( nodes[ j-1 ]->GetPosition()->GetTypeOfPosition() == SMDS_TOP_3DSPACE )
-          {
-            const PointGeomInfo& pgi = elem.GeomInfoPi(j);
-            meshDS->SetNodeOnFace( nodes[ j-1 ], faceID, pgi.u, pgi.v);
-          }
-        }
-        if ( j > elem.GetNP() )
-        {
-          if ( elem.GetType() == TRIG )
-            helper.AddFace(nodes[0],nodes[1],nodes[2]);
-          else
-            helper.AddFace(nodes[0],nodes[1],nodes[2],nodes[3]);
-        }
-      }
+      FillNodesAndElements( aMesh, helper, ngMesh, nodeVec, faceID );      
 
       break;
     } // two attempts
index 7d103ed4e7a6d1fbddf708435089744734ea9273..1be47c526850b8e61d6ec35b95292721512b4f6e 100644 (file)
 
 #include <SMESH_Algo.hxx>
 #include <SMESH_Mesh.hxx>
+#include <SMESH_Group.hxx>
+#include <SMESHDS_GroupBase.hxx>
+
+#include "NETGENPlugin_Mesher.hxx"
 
 class StdMeshers_MaxElementArea;
 class StdMeshers_LengthFromEdges;
@@ -50,7 +54,7 @@ public:
                                Hypothesis_Status&  aStatus);
 
   virtual bool Compute(SMESH_Mesh&         aMesh,
-                       const TopoDS_Shape& aShape);
+                       const TopoDS_Shape& aShape);  
 
   virtual void CancelCompute();
 
@@ -59,6 +63,29 @@ public:
   virtual bool Evaluate(SMESH_Mesh& aMesh, const TopoDS_Shape& aShape,
                         MapShapeNbElems& aResMap);
 
+  bool MapSegmentsToEdges(SMESH_Mesh&         aMesh,
+                          const TopoDS_Shape& aShape, 
+                          NETGENPlugin_NetgenLibWrapper &ngLib, 
+                          vector< const SMDS_MeshNode* >& nodeVec, 
+                          std::map<int,const SMDS_MeshNode*>& premeshedNodes, 
+                          std::map<int,std::vector<double>>& newNetgenCoordinates,
+                          std::map<int,std::vector<smIdType>>& newNetgenElements );  
+
+  std::tuple<bool,bool> SetParameteres( SMESH_Mesh& aMesh, const TopoDS_Shape& aShape, 
+                                        NETGENPlugin_Mesher& aMesher,  netgen::Mesh * ngMeshes,
+                                        netgen::OCCGeometry& occgeoComm, bool isSubMeshSupported = true );     
+
+  bool ComputeMaxhOfFace( TopoDS_Face& Face, NETGENPlugin_Mesher& aMesher, TSideVector& wires, 
+                          netgen::OCCGeometry& occgeoComm, bool isDefaultHyp, bool isCommonLocalSize );
+  
+  void FillNodesAndElements( SMESH_Mesh& aMesh, SMESH_MesherHelper& helper, netgen::Mesh * ngMesh, vector< const SMDS_MeshNode* >& nodeVec, int faceId );
+
+ /*!
+ * \brief FillNodesAndElements, fill created triangular elements by netgen to the smesh data structure
+ */
+  void FillNodesAndElements( SMESH_Mesh& aMesh, SMESH_MesherHelper& helper, netgen::Mesh * ngMesh, vector< const SMDS_MeshNode* >& nodeVec, map<int, const SMDS_MeshNode* >& ng2smesh,
+                              std::map<int,std::vector<double>>& newNetgenCoordinates, std::map<int,std::vector<smIdType>>& newNetgenElements, const int numberOfPremeshedNodes );
+
 protected:
   const StdMeshers_MaxElementArea*       _hypMaxElementArea;
   const StdMeshers_LengthFromEdges*      _hypLengthFromEdges;
diff --git a/src/NETGENPlugin/NETGENPlugin_NETGEN_2D_Remote.cxx b/src/NETGENPlugin/NETGENPlugin_NETGEN_2D_Remote.cxx
new file mode 100644 (file)
index 0000000..5e870e8
--- /dev/null
@@ -0,0 +1,352 @@
+// Copyright (C) 2007-2023  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
+//
+
+//=============================================================================
+// File      : NETGENPlugin_NETGEN_2D_Remote.cxx
+// Created   : mardi 12 Decembre 2023
+// Author    : Cesar Conopoima (OCC)
+// Project   : SALOME
+//=============================================================================
+//
+//
+#include "NETGENPlugin_NETGEN_2D_Remote.hxx"
+
+#include "NETGENPlugin_DriverParam.hxx"
+#include "NETGENPlugin_Hypothesis.hxx"
+
+#include "Utils_SALOME_Exception.hxx"
+
+#include <SMESH_Gen.hxx>
+#include <SMESH_Mesh.hxx>
+#include <SMESH_ParallelMesh.hxx>
+#include <SMESH_MesherHelper.hxx>
+#include <SMESH_DriverShape.hxx>
+#include <SMESH_DriverMesh.hxx>
+#include <SMESHDS_Mesh.hxx>
+#include <SMESH_MeshLocker.hxx>
+
+#include <QString>
+#include <QProcess>
+
+#include <boost/filesystem.hpp>
+namespace fs = boost::filesystem;
+/*
+  Netgen include files
+*/
+
+#ifndef OCCGEOMETRY
+#define OCCGEOMETRY
+#endif
+#include <occgeom.hpp>
+
+#ifdef NETGEN_V5
+#include <ngexception.hpp>
+#endif
+#ifdef NETGEN_V6
+#include <core/exception.hpp>
+#endif
+
+namespace nglib {
+#include <nglib.h>
+}
+namespace netgen {
+
+  NETGENPLUGIN_DLL_HEADER
+  extern MeshingParameters mparam;
+
+  NETGENPLUGIN_DLL_HEADER
+  extern volatile multithreadt multithread;
+}
+using namespace nglib;
+
+//=============================================================================
+/*!
+ * Constructor
+ */
+//=============================================================================
+
+NETGENPlugin_NETGEN_2D_Remote::NETGENPlugin_NETGEN_2D_Remote(int hypId, SMESH_Gen * gen)
+  : NETGENPlugin_NETGEN_2D_ONLY(hypId, gen)
+{
+  _name = "NETGEN_2D_Remote";
+}
+
+//=============================================================================
+/*!
+ * Destructor
+ */
+//=============================================================================
+
+NETGENPlugin_NETGEN_2D_Remote::~NETGENPlugin_NETGEN_2D_Remote()
+{
+}
+
+/**
+ * @brief Fill the structure netgen_param with the information from the hypothesis
+ *
+ * @param hyp the hypothesis
+ * @param aParams the netgen_param structure
+ */
+void NETGENPlugin_NETGEN_2D_Remote::fillParameters(const NETGENPlugin_Hypothesis* hyp, netgen_params &aParams)
+{
+  aParams.myType             = hypoType::Hypo;
+  aParams.maxh               = hyp->GetMaxSize();
+  aParams.minh               = hyp->GetMinSize();
+  aParams.segmentsperedge    = hyp->GetNbSegPerEdge();
+  aParams.grading            = hyp->GetGrowthRate();
+  aParams.curvaturesafety    = hyp->GetNbSegPerRadius();
+  aParams.secondorder        = hyp->GetSecondOrder() ? 1 : 0;
+  aParams.quad               = hyp->GetQuadAllowed() ? 1 : 0;
+  aParams.optimize           = hyp->GetOptimize();
+  aParams.fineness           = hyp->GetFineness();
+  aParams.uselocalh          = hyp->GetSurfaceCurvature();
+  aParams.merge_solids       = hyp->GetFuseEdges();
+  aParams.chordalError       = hyp->GetChordalErrorEnabled() ? hyp->GetChordalError() : -1.;
+  aParams.optsteps2d         = aParams.optimize ? hyp->GetNbSurfOptSteps() : 0;
+  aParams.optsteps3d         = aParams.optimize ? hyp->GetNbVolOptSteps()  : 0;
+  aParams.elsizeweight       = hyp->GetElemSizeWeight();
+  aParams.opterrpow          = hyp->GetWorstElemMeasure();
+  aParams.delaunay           = hyp->GetUseDelauney();
+  aParams.checkoverlap       = hyp->GetCheckOverlapping();
+  aParams.checkchartboundary = hyp->GetCheckChartBoundary();
+#ifdef NETGEN_V6
+  // std::string
+  aParams.meshsizefilename = hyp->GetMeshSizeFile();
+  aParams.closeedgefac = 2;
+  aParams.nbThreads = hyp->GetNbThreads();
+#else
+  // const char*
+  aParams.meshsizefilename = hyp->GetMeshSizeFile();
+  aParams.closeedgefac = 0;
+  aParams.nbThreads = 0;
+#endif
+}
+
+/**
+ * @brief write in a binary file the orientation for each surface element of the mesh
+ *
+ * @param aMesh The mesh
+ * @param aShape the shape associated to the mesh
+ * @param output_file name of the binary file
+ */
+void NETGENPlugin_NETGEN_2D_Remote::exportElementOrientation(SMESH_Mesh& aMesh,
+                                                      const TopoDS_Shape& aShape,
+                                                      const std::string output_file)
+{
+  std::ofstream df(output_file, ios::out|ios::binary);
+  int size=0;
+  df.write((char*)&size, sizeof(int));
+  df.close();
+}
+
+/**
+ * @brief Compute mesh associate to shape
+ *
+ * @param aMesh The mesh
+ * @param aShape The shape
+ * @return true fi there are some error
+ */
+bool NETGENPlugin_NETGEN_2D_Remote::Compute(SMESH_Mesh&         aMesh,
+                                           const TopoDS_Shape& aShape)
+{
+  {
+    SMESH_MeshLocker myLocker(&aMesh);
+    SMESH_Hypothesis::Hypothesis_Status hypStatus;
+    NETGENPlugin_NETGEN_2D_ONLY::CheckHypothesis(aMesh, aShape, hypStatus);
+  }
+  SMESH_ParallelMesh& aParMesh = dynamic_cast<SMESH_ParallelMesh&>(aMesh);
+
+  // Temporary folder for run
+#ifdef WIN32
+  fs::path tmp_folder = aParMesh.GetTmpFolder() / fs::path("Volume-%%%%-%%%%");
+#else
+  fs::path tmp_folder = aParMesh.GetTmpFolder() / fs::unique_path(fs::path("Volume-%%%%-%%%%"));
+#endif
+  fs::create_directories(tmp_folder);
+  // Using MESH2D generated after all triangles where created.
+  fs::path mesh_file=aParMesh.GetTmpFolder() / fs::path("Mesh1D.med"); // read the premeshed elements from 2D version
+  fs::path element_orientation_file=tmp_folder / fs::path("element_orientation.dat");
+  fs::path new_element_file=tmp_folder / fs::path("new_elements.dat");
+  // Not used kept for debug
+  //fs::path output_mesh_file=tmp_folder / fs::path("output_mesh.med");
+  fs::path shape_file=tmp_folder / fs::path("shape.brep");
+  fs::path param_file=tmp_folder / fs::path("netgen_lenghtfromedge.txt"); /*becuase name contain 'lenghtfromedge' set length of 2D from premeshed 1D elements*/
+  fs::path log_file=tmp_folder / fs::path("run.log");
+  fs::path cmd_file=tmp_folder / fs::path("cmd.txt");
+  std::string mesh_name = "MESH";
+
+  {
+    SMESH_MeshLocker myLocker(&aMesh);
+    //Writing Shape
+    SMESH_DriverShape::exportShape(shape_file.string(), aShape);
+
+    //Writing hypo
+    // netgen_params aParams;
+    // fillParameters(_hypParameters, aParams);
+    // exportNetgenParams(param_file.string(), aParams);
+    {
+      // Simply write the file with the proper name
+      std::ofstream myfile(param_file);
+      myfile << 1 << std::endl;
+      myfile.close();
+    }
+      
+    // Exporting element orientation
+    exportElementOrientation(aMesh, aShape, element_orientation_file.string());
+  }
+
+  // Calling run_mesher
+  // Path to mesher script
+  fs::path mesher_launcher = fs::path(std::getenv("SMESH_ROOT_DIR"))/
+       fs::path("bin")/
+       fs::path("salome")/
+       fs::path("mesher_launcher.py");
+
+
+  std::string s_program="python3";
+  std::list<std::string> params;
+  params.push_back(mesher_launcher.string());
+  params.push_back("NETGEN2D");
+  params.push_back(mesh_file.string());
+  params.push_back(shape_file.string());
+  params.push_back(param_file.string());
+  params.push_back("--elem-orient-file=" + element_orientation_file.string());
+  params.push_back("--new-element-file=" + new_element_file.string());
+  // params.push_back("--output-mesh-file=" + output_mesh_file.string());
+
+  // Parallelism method parameters
+  int method = aParMesh.GetParallelismMethod();
+  if(method == ParallelismMethod::MultiThread){
+    params.push_back("--method=local");
+  } else if (method == ParallelismMethod::MultiNode){
+    params.push_back("--method=cluster");
+    params.push_back("--resource="+aParMesh.GetResource());
+    params.push_back("--wc-key="+aParMesh.GetWcKey());
+    params.push_back("--nb-proc=1");
+    params.push_back("--nb-proc-per-node="+to_string(aParMesh.GetNbProcPerNode()));
+    params.push_back("--nb-node="+to_string(aParMesh.GetNbNode()));
+    params.push_back("--walltime="+aParMesh.GetWalltime());
+  } else {
+    throw SALOME_Exception("Unknown parallelism method "+method);
+  }
+  std::string cmd = "";
+  cmd += s_program;
+  for(auto arg: params){
+    cmd += " " + arg;
+  }
+  MESSAGE("Running command: ");
+  MESSAGE(cmd);
+  // Writing command in cmd.log
+  {
+    std::ofstream flog(cmd_file.string());
+    flog << cmd << endl;
+  }
+
+  // Building arguments for QProcess
+  QString program = QString::fromStdString(s_program);
+  QStringList arguments;
+  for(auto arg : params){
+    arguments << arg.c_str();
+  }
+
+  QString out_file = log_file.string().c_str();
+  QProcess myProcess;
+  // myProcess.setProcessChannelMode(QProcess::MergedChannels);
+  myProcess.setProcessChannelMode(QProcess::ForwardedChannels);
+  myProcess.setStandardOutputFile(out_file);
+
+  myProcess.start(program, arguments);
+  // Waiting for process to finish (argument -1 make it wait until the end of
+  // the process otherwise it just waits 30 seconds)
+  bool finished = myProcess.waitForFinished(-1);
+  int ret = myProcess.exitCode();
+  if(ret != 0 || !finished){
+    // Run crahed
+    std::string msg = "Issue with mesh_launcher: \n";
+    msg += "See log for more details: " + log_file.string() + "\n";
+    msg += cmd + "\n";
+    throw SALOME_Exception(msg);
+  }
+
+  {
+    SMESH_MeshLocker myLocker(&aMesh);
+    std::ifstream df(new_element_file.string(), ios::binary);
+
+    int totalPremeshedNodes;
+    int NetgenNbOfNodes;
+    int NetgenNbOfNodesNew;
+    int NetgenNbOfTriangles;
+    double NetgenPoint[3];
+    int    NetgenTriangle[3];
+    int nodeID;
+
+    SMESH_MesherHelper helper(aMesh);
+    // This function is mandatory for setElementsOnShape to work
+    helper.IsQuadraticSubMesh(aShape);
+    helper.SetElementsOnShape( true );
+
+    df.read((char*) &totalPremeshedNodes, sizeof(int));
+    // Number of nodes in intial mesh
+    df.read((char*) &NetgenNbOfNodes, sizeof(int));
+    // Number of nodes added by netgen
+    df.read((char*) &NetgenNbOfNodesNew, sizeof(int));
+
+    // Filling nodevec (correspondence netgen numbering mesh numbering)
+    vector< const SMDS_MeshNode* > nodeVec ( NetgenNbOfNodesNew + 2 );
+    SMESHDS_Mesh * meshDS = helper.GetMeshDS();
+    for (int nodeIndex = 1; nodeIndex <= NetgenNbOfNodes; ++nodeIndex )
+    {
+      //Id of the point
+      df.read((char*) &nodeID, sizeof(int));
+      nodeVec.at(nodeID) = meshDS->FindNode(nodeID);
+    }
+
+    // Add new points and update nodeVec
+    for (int nodeIndex = totalPremeshedNodes + 1; nodeIndex <= NetgenNbOfNodesNew; ++nodeIndex )
+    {
+      df.read((char *) &NetgenPoint, sizeof(double)*3);
+      nodeVec.at(nodeIndex) = helper.AddNode(NetgenPoint[0], NetgenPoint[1], NetgenPoint[2]);
+    }
+
+    // Add triangles
+    df.read((char*) &NetgenNbOfTriangles, sizeof(int));
+    for ( int elemIndex = 1; elemIndex <= NetgenNbOfTriangles; ++elemIndex )
+    {
+      df.read((char*) &NetgenTriangle, sizeof(int)*3);
+      if ( nodeVec.at( NetgenTriangle[0] ) && nodeVec.at( NetgenTriangle[1] ) && nodeVec.at( NetgenTriangle[2] ) )
+        helper.AddFace(nodeVec.at( NetgenTriangle[0] ), nodeVec.at( NetgenTriangle[1] ), nodeVec.at( NetgenTriangle[2] ) );            
+    }
+  }
+
+  return true;
+}
+
+/**
+ * @brief Assign submeshes to compute
+ *
+ * @param aSubMesh submesh to add
+ */
+void NETGENPlugin_NETGEN_2D_Remote::setSubMeshesToCompute(SMESH_subMesh * aSubMesh)
+{
+  SMESH_MeshLocker myLocker(aSubMesh->GetFather());
+  SMESH_Algo::setSubMeshesToCompute(aSubMesh);
+}
diff --git a/src/NETGENPlugin/NETGENPlugin_NETGEN_2D_Remote.hxx b/src/NETGENPlugin/NETGENPlugin_NETGEN_2D_Remote.hxx
new file mode 100644 (file)
index 0000000..4e90666
--- /dev/null
@@ -0,0 +1,75 @@
+// Copyright (C) 2007-2023  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
+//
+
+//=============================================================================
+// File      : NETGENPlugin_NETGEN_2D_Remote.hxx
+// Created   : mardi 12 Decembre 2023
+// Author    : Cesar Conopoima (OCC)
+// Project   : SALOME
+//=============================================================================
+//
+#ifndef _NETGENPlugin_NETGEN_2D_REMOTE_HXX_
+#define _NETGENPlugin_NETGEN_2D_REMOTE_HXX_
+
+#include "NETGENPlugin_NETGEN_2D_ONLY.hxx"
+
+#include <vector>
+#include <map>
+
+class StdMeshers_ViscousLayers;
+class StdMeshers_MaxElementVolume;
+class NETGENPlugin_Hypothesis;
+class NETGENPlugin_NetgenLibWrapper;
+class netgen_params;
+class SMDS_MeshNode;
+
+using namespace std;
+
+class NETGENPLUGIN_EXPORT NETGENPlugin_NETGEN_2D_Remote: public NETGENPlugin_NETGEN_2D_ONLY
+{
+ public:
+  NETGENPlugin_NETGEN_2D_Remote(int hypId, SMESH_Gen* gen);
+  virtual ~NETGENPlugin_NETGEN_2D_Remote();
+
+  // Function whould not be used with remote Computing
+  bool CheckHypothesis (SMESH_Mesh&         aMesh,
+                        const TopoDS_Shape& aShape,
+                        Hypothesis_Status&  aStatus) override {(void)aMesh;(void)aShape;aStatus = HYP_OK;return true;};
+
+  bool Compute(SMESH_Mesh&         aMesh,
+              const TopoDS_Shape& aShape) override;
+
+  void setSubMeshesToCompute(SMESH_subMesh * aSubMesh) override;
+
+
+ protected:
+ void exportElementOrientation(SMESH_Mesh& aMesh,
+                                const TopoDS_Shape& aShape,
+                                const std::string output_file);
+                                
+  void fillParameters(const NETGENPlugin_Hypothesis* hyp,
+                      netgen_params &aParams);
+
+
+};
+
+#endif
diff --git a/src/NETGENPlugin/NETGENPlugin_NETGEN_2D_Remote_i.cxx b/src/NETGENPlugin/NETGENPlugin_NETGEN_2D_Remote_i.cxx
new file mode 100644 (file)
index 0000000..3aee667
--- /dev/null
@@ -0,0 +1,77 @@
+// Copyright (C) 2007-2023  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
+//
+
+//  SMESH SMESH_I : idl implementation based on 'SMESH' unit's classes
+//  File   : NETGENPlugin_NETGEN_2D_Remote_i.cxx
+//  Author : Cesar Conopoima (OCC)
+//  Module : NETGENPlugin
+//  $Header$
+//
+#include "NETGENPlugin_NETGEN_2D_Remote_i.hxx"
+#include "SMESH_Gen.hxx"
+
+#include "Utils_CorbaException.hxx"
+#include "utilities.h"
+
+//=============================================================================
+/*!
+ *  NETGENPlugin_NETGEN_2D_Remote_i::NETGENPlugin_NETGEN_2D_Remote_i
+ *
+ *  Constructor
+ */
+//=============================================================================
+
+NETGENPlugin_NETGEN_2D_Remote_i::NETGENPlugin_NETGEN_2D_Remote_i( PortableServer::POA_ptr thePOA,
+                                                    ::SMESH_Gen*            theGenImpl )
+     : SALOME::GenericObj_i( thePOA ),
+       SMESH_Hypothesis_i( thePOA ),
+       SMESH_Algo_i( thePOA ),
+       SMESH_2D_Algo_i( thePOA )
+{
+  myBaseImpl = new ::NETGENPlugin_NETGEN_2D_Remote( theGenImpl->GetANewId(), theGenImpl );
+}
+
+//=============================================================================
+/*!
+ *  NETGENPlugin_NETGEN_2D_Remote_i::~NETGENPlugin_NETGEN_2D_Remote_i
+ *
+ *  Destructor
+ */
+//=============================================================================
+
+NETGENPlugin_NETGEN_2D_Remote_i::~NETGENPlugin_NETGEN_2D_Remote_i()
+{
+}
+
+//=============================================================================
+/*!
+ *  NETGENPlugin_NETGEN_2D_Remote_i::GetImpl
+ *
+ *  Get implementation
+ */
+//=============================================================================
+
+::NETGENPlugin_NETGEN_2D_Remote* NETGENPlugin_NETGEN_2D_Remote_i::GetImpl()
+{
+  return ( ::NETGENPlugin_NETGEN_2D_Remote* )myBaseImpl;
+}
+
diff --git a/src/NETGENPlugin/NETGENPlugin_NETGEN_2D_Remote_i.hxx b/src/NETGENPlugin/NETGENPlugin_NETGEN_2D_Remote_i.hxx
new file mode 100644 (file)
index 0000000..654637f
--- /dev/null
@@ -0,0 +1,58 @@
+// Copyright (C) 2007-2023  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
+//
+
+//  SMESH SMESH_I : idl implementation based on 'SMESH' unit's classes
+//  File   : NETGENPlugin_NETGEN_2D_Remote_i.hxx
+//  Author : Cesar Conopoima (OCC)
+//  Module : NETGENPlugin
+//  $Header$
+//
+#ifndef _NETGENPlugin_NETGEN_2D_REMOTE_I_HXX_
+#define _NETGENPlugin_NETGEN_2D_REMOTE_I_HXX_
+
+#include "NETGENPlugin_Defs.hxx"
+
+#include <SALOMEconfig.h>
+#include CORBA_SERVER_HEADER(NETGENPlugin_Algorithm)
+
+#include "SMESH_2D_Algo_i.hxx"
+#include "NETGENPlugin_NETGEN_2D_Remote.hxx"
+
+// ======================================================
+// NETGEN 2d algorithm
+// ======================================================
+class NETGENPLUGIN_EXPORT NETGENPlugin_NETGEN_2D_Remote_i:
+  public virtual POA_NETGENPlugin::NETGENPlugin_NETGEN_2D_Remote,
+  public virtual SMESH_2D_Algo_i
+{
+public:
+  // Constructor
+  NETGENPlugin_NETGEN_2D_Remote_i( PortableServer::POA_ptr thePOA,
+                            ::SMESH_Gen*            theGenImpl );
+  // Destructor
+  virtual ~NETGENPlugin_NETGEN_2D_Remote_i();
+
+  // Get implementation
+  ::NETGENPlugin_NETGEN_2D_Remote* GetImpl();
+};
+
+#endif
diff --git a/src/NETGENPlugin/NETGENPlugin_NETGEN_2D_SA.cxx b/src/NETGENPlugin/NETGENPlugin_NETGEN_2D_SA.cxx
new file mode 100644 (file)
index 0000000..d06d4e9
--- /dev/null
@@ -0,0 +1,270 @@
+// Copyright (C) 2007-2023  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
+//
+
+// NETGENPlugin : C++ implementation
+// File   : NETGENPlugin_NETGEN_2D_SA.cxx
+// Author : Cesar Conopoima (OCC)
+// Date   : 23/10/2023
+// Project   : SALOME
+//=============================================================================
+//
+#include "NETGENPlugin_NETGEN_2D_SA.hxx"
+#include "NETGENPlugin_Hypothesis_2D.hxx"
+#include "NETGENPlugin_SimpleHypothesis_2D.hxx"
+#include "NETGENPlugin_Mesher.hxx"
+
+#include <SMESHDS_Mesh.hxx>
+#include <SMESH_ControlsDef.hxx>
+#include <SMESH_Gen.hxx>
+#include <SMESH_Mesh.hxx>
+#include <SMESH_MesherHelper.hxx>
+#include <SMESH_DriverShape.hxx>
+#include <SMESH_DriverMesh.hxx>
+#include <StdMeshers_LengthFromEdges.hxx>
+#include <StdMeshers_MaxElementArea.hxx>
+#include <utilities.h>
+
+#ifdef WIN32
+#include <filesystem>
+namespace fs = std::filesystem;
+#else
+#include <boost/filesystem.hpp>
+namespace fs = boost::filesystem;
+#endif
+
+namespace nglib {
+#include <nglib.h>
+}
+#include <meshing.hpp>
+
+using namespace nglib;
+
+//=============================================================================
+/*!
+ *  
+ */
+//=============================================================================
+
+NETGENPlugin_NETGEN_2D_SA::NETGENPlugin_NETGEN_2D_SA()
+  : NETGENPlugin_NETGEN_2D_ONLY(0, new SMESH_Gen())
+{
+  _name = "NETGEN_2D_SA";  
+}
+
+//=============================================================================
+/*!
+ *  
+ */
+//=============================================================================
+
+NETGENPlugin_NETGEN_2D_SA::~NETGENPlugin_NETGEN_2D_SA()
+{
+}
+
+/**
+ * @brief fill plugin hypothesis from the netgen_params structure
+ *
+ * @param aParams the structure
+ */
+void NETGENPlugin_NETGEN_2D_SA::fillHyp(const std::string param_file, netgen_params aParams)
+{  
+
+  if( aParams.has_netgen_param && aParams.myType == hypoType::Hypo )
+  {
+    NETGENPlugin_Hypothesis_2D * hyp = new NETGENPlugin_Hypothesis_2D(0, GetGen());
+    hyp->SetMaxSize(aParams.maxh);
+    hyp->SetMinSize(aParams.minh);
+    hyp->SetGrowthRate(aParams.grading);
+    hyp->SetQuadAllowed(aParams.quad);
+    hyp->SetFineness((NETGENPlugin_Hypothesis::Fineness)aParams.fineness);
+    hyp->SetChordalErrorEnabled(aParams.chordalError);
+    if(aParams.optimize){
+      hyp->SetNbSurfOptSteps(aParams.optsteps2d);
+    }
+    _hypParameters = dynamic_cast< const NETGENPlugin_Hypothesis_2D *> (hyp);
+  }
+  else 
+  {
+    NETGENPlugin_Hypothesis * hyp = new NETGENPlugin_Hypothesis(0, GetGen());
+    if (param_file.find("lenghtfromedge") != std::string::npos )
+    {
+      _hypLengthFromEdges = dynamic_cast<const StdMeshers_LengthFromEdges*> (hyp);
+    }
+    else if ( param_file.find("maxarea") != std::string::npos )
+    {
+      StdMeshers_MaxElementArea * hyp = new StdMeshers_MaxElementArea(0, GetGen());
+
+      std::ifstream myfile(param_file);
+      std::string line;
+      double maxArea;
+      std::getline(myfile, line);
+      bool hashypothesis = std::stoi(line);
+      if ( hashypothesis )
+      {
+        std::getline(myfile, line);
+        maxArea = std::stod(line);
+        hyp->SetMaxArea( maxArea );
+        _hypMaxElementArea = static_cast<const StdMeshers_MaxElementArea*> (hyp);
+      }
+      
+      myfile.close();
+    }
+  }  
+}
+
+
+/**
+ * @brief Write a binary file containing information on the elements/nodes
+ *        created by the mesher
+ *
+ * @param premeshedNodes map of the premeshed nodes of 1D elements
+ * @param newNetgenCoordinates map of the coordinate of new netgen points created in the mesh
+ * @param newNetgenElements map of the element triangulation
+ * @return true if there are some error
+ */
+bool NETGENPlugin_NETGEN_2D_SA::fillNewElementFile( std::string new_element_file,
+                                                    const int numberOfGlobalPremeshedNodes,
+                                                    std::map<int,const SMDS_MeshNode*>& premeshedNodes, 
+                                                    std::map<int,std::vector<double>>& newNetgenCoordinates,
+                                                    std::map<int,std::vector<smIdType>>& newNetgenElements )
+{
+  MESSAGE("Writting new elements")
+
+  int NetgenNbOfNodes     = premeshedNodes.size();
+  int NetgenNbOfNodesNew  = numberOfGlobalPremeshedNodes + newNetgenCoordinates.size();
+  int NetgenNbOfTriangles = newNetgenElements.size();
+  bool isOK = ( NetgenNbOfTriangles > 0 );
+  if ( isOK && !new_element_file.empty() )
+  {
+    int    NetgenElement[3];
+    std::ofstream df(new_element_file, ios::out|ios::binary);
+    // Writing nodevec (correspondence netgen numbering mesh numbering)
+    // Number of nodes
+    df.write((char*) &numberOfGlobalPremeshedNodes, sizeof(int));
+    df.write((char*) &NetgenNbOfNodes, sizeof(int));
+    df.write((char*) &NetgenNbOfNodesNew, sizeof(int));
+    for (auto k : premeshedNodes )
+      df.write((char*) &k.first, sizeof(int));      
+        
+    // Writing info on new points
+    for (auto k : newNetgenCoordinates )
+      df.write((char*) k.second.data(), sizeof(double)*3 ); 
+      
+    // create triangles
+    df.write((char*) &NetgenNbOfTriangles, sizeof(int));
+    for ( int elemIndex = 1; elemIndex <= NetgenNbOfTriangles; ++elemIndex )
+    {
+      int nodes = newNetgenElements[ elemIndex ].size();
+      for (int i = 0; i < nodes; i++)
+        NetgenElement[ i ] = (int) newNetgenElements[ elemIndex ][ i ];
+      
+      df.write((char*) &NetgenElement, sizeof(int)* nodes );
+    }
+
+    df.close();
+  }
+  return false;
+}
+
+/**
+ * @brief Compute the mesh based on the 
+ *
+ * @param input_mesh_file Mesh file (containing 2D elements)
+ * @param shape_file Shape file (BREP or STEP format)
+ * @param hypo_file Ascii file containing the netgen parameters
+ * @param element_orientation_file Binary file containing the orientation of surface elemnts
+ * @param new_element_file output file containing info the elements created by the mesher
+ * @param output_mesh_file output mesh file (if empty it will not be created)
+ * @return true, false
+ */
+bool NETGENPlugin_NETGEN_2D_SA::Compute( SMESH_Mesh& aMesh, TopoDS_Shape &aShape, std::string new_element_file )
+{
+  // Nodes on edge are double because each face is treat one by one, so seam edges are iterated twice 
+  vector< const SMDS_MeshNode* > nodeVec;
+  NETGENPlugin_Mesher aMesher( &aMesh, aShape, false /*isVolume=*/ );
+  NETGENPlugin_NetgenLibWrapper ngLib;
+  ngLib._isComputeOk = false;
+
+  fs::path netgen_log_file = fs::path(new_element_file).remove_filename() / fs::path("NETGEN.out");
+  MESSAGE("netgen ouput"<<netgen_log_file.string());
+  ngLib.setOutputFile(netgen_log_file.string());
+  
+  netgen::OCCGeometry occgeoComm;
+  NETGENPlugin_NETGEN_2D_ONLY::SetParameteres( aMesh, aShape, aMesher, ngLib._ngMesh, occgeoComm, false /*submesh is not supported*/ );
+
+  std::map<int,const SMDS_MeshNode*> premeshedNodes;
+  std::map<int,std::vector<double>> newNetgenCoordinates;
+  std::map<int,std::vector<smIdType>> newNetgenElements;
+  const int numberOfTotalPremeshedNodes = aMesh.NbNodes();
+  bool compute = NETGENPlugin_NETGEN_2D_ONLY::MapSegmentsToEdges( aMesh, aShape, ngLib, nodeVec,
+                                                                  premeshedNodes, newNetgenCoordinates, 
+                                                                  newNetgenElements );
+  
+  compute = fillNewElementFile(new_element_file, 
+                                numberOfTotalPremeshedNodes, 
+                                premeshedNodes, 
+                                newNetgenCoordinates, 
+                                newNetgenElements);
+  return compute;
+}
+
+/**
+ * @brief Running the mesher on the given files
+ *
+ * @param input_mesh_file Mesh file (containing 2D elements)
+ * @param shape_file Shape file (BREP or STEP format)
+ * @param hypo_file Ascii file containing the netgen parameters
+ * @param element_orientation_file Binary file containing the orientation of surface elemnts
+ * @param new_element_file output file containing info the elements created by the mesher
+ * @param output_mesh_file output mesh file (if empty it will not be created)
+ * @return int
+ */
+int NETGENPlugin_NETGEN_2D_SA::run(const std::string input_mesh_file,
+                                    const std::string shape_file,
+                                    const std::string hypo_file,
+                                    const std::string element_orientation_file,
+                                    const std::string new_element_file,
+                                    const std::string output_mesh_file)
+{
+
+  std::unique_ptr<SMESH_Mesh> myMesh(_gen->CreateMesh(false));
+
+  SMESH_DriverMesh::importMesh(input_mesh_file, *myMesh);
+  // Importing shape
+  TopoDS_Shape myShape;
+  SMESH_DriverShape::importShape(shape_file, myShape);
+  // Importing hypothesis
+  netgen_params myParams;  
+  importNetgenParams(hypo_file, myParams);
+  fillHyp(hypo_file,myParams);  
+  MESSAGE("Meshing with netgen2d");
+  int ret = (int) Compute( *myMesh, myShape, new_element_file );
+
+  if(ret){
+    std::cerr << "Meshing failed" << std::endl;
+    return ret;
+  }
+
+  if(!output_mesh_file.empty()){
+    std::string meshName = "MESH";
+    SMESH_DriverMesh::exportMesh(output_mesh_file, *myMesh, meshName);
+  }
+
+  return ret;
+}
diff --git a/src/NETGENPlugin/NETGENPlugin_NETGEN_2D_SA.hxx b/src/NETGENPlugin/NETGENPlugin_NETGEN_2D_SA.hxx
new file mode 100644 (file)
index 0000000..99e4706
--- /dev/null
@@ -0,0 +1,63 @@
+// Copyright (C) 2007-2023  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
+//
+
+// NETGENPlugin : C++ implementation
+// File      : NETGENPlugin_NETGEN_2D_SA.hxx
+// Author    : Cesar Conopoima (OCC)
+// Date      : 23/10/2023
+// Project   : SALOME
+//=============================================================================
+//
+#ifndef _NETGENPlugin_NETGEN_2D_SA_HXX_
+#define _NETGENPlugin_NETGEN_2D_SA_HXX_
+
+#include "NETGENPlugin_Defs.hxx"
+#include "NETGENPlugin_DriverParam.hxx"
+#include "NETGENPlugin_NETGEN_2D_ONLY.hxx"
+#include "SMESH_Algo.hxx"
+#include "SMESH_Mesh.hxx"
+
+class NETGENPlugin_Mesher;
+
+class NETGENPLUGIN_EXPORT NETGENPlugin_NETGEN_2D_SA: public NETGENPlugin_NETGEN_2D_ONLY
+{
+public:
+  NETGENPlugin_NETGEN_2D_SA();
+  virtual ~NETGENPlugin_NETGEN_2D_SA();
+
+  bool Compute(SMESH_Mesh& aMesh, TopoDS_Shape &aShape, std::string new_element_file );
+  
+  int run(const std::string input_mesh_file,
+          const std::string shape_file,
+          const std::string hypo_file,
+          const std::string element_orientation_file,
+          const std::string new_element_file,
+          const std::string output_mesh_file);
+
+  bool fillNewElementFile( std::string new_element_file, 
+                           const int numberOfGlobalPremeshedNodes,
+                           std::map<int,const SMDS_MeshNode*>& premeshedNodes, 
+                           std::map<int,std::vector<double>>& newNetgenCoordinates,
+                           std::map<int,std::vector<smIdType>>& newNetgenElements );
+protected:
+  
+  void fillHyp(const std::string param_file, netgen_params aParams);
+};
+
+#endif
index d5eb93eb06261c5475c78cca83056a14bff4ee98..7f441983e71e1f58c4d37d5983b1f53edabd203f 100644 (file)
@@ -437,7 +437,6 @@ bool NETGENPlugin_NETGEN_3D::computePrepareParam(
 
   NETGENPlugin_Mesher aMesher( &aMesh, helper.GetSubShape(), /*isVolume=*/true );
 
-
   if ( _hypParameters )
   {
     aMesher.SetParameters( _hypParameters );
@@ -637,7 +636,6 @@ bool NETGENPlugin_NETGEN_3D::Compute(
   computeRunMesher(occgeo, nodeVec, ngLib._ngMesh, ngLib, startWith, endWith);
 
   computeFillMesh(nodeVec, ngLib, helper, Netgen_NbOfNodes);
-
   return false;
 
 }
index 6dcda6fed74306ebcb27ca0bceaec04d071a6532..a84d16607bd3ab66566ec5afd81c6256590915da 100644 (file)
@@ -110,6 +110,7 @@ NETGENPlugin_NETGEN_3D_Remote::~NETGENPlugin_NETGEN_3D_Remote()
  */
 void NETGENPlugin_NETGEN_3D_Remote::fillParameters(const NETGENPlugin_Hypothesis* hyp, netgen_params &aParams)
 {
+  aParams.myType             = hypoType::Hypo;
   aParams.maxh               = hyp->GetMaxSize();
   aParams.minh               = hyp->GetMinSize();
   aParams.segmentsperedge    = hyp->GetNbSegPerEdge();
index 75ee83431c739331da20992b7128afa8ed318329..135eb0a08b9bcda2fe85a43e98a74d71aac66b6d 100644 (file)
@@ -168,7 +168,6 @@ bool NETGENPlugin_NETGEN_3D_SA::computeFillNewElementFile(
 
   int Netgen_NbOfNodesNew = Ng_GetNP(Netgen_mesh);
   int Netgen_NbOfTetra    = Ng_GetNE(Netgen_mesh);
-
   bool isOK = ( Netgen_NbOfTetra > 0 );
   if ( isOK && !new_element_file.empty() )
   {
@@ -289,7 +288,6 @@ int NETGENPlugin_NETGEN_3D_SA::run(const std::string input_mesh_file,
 
   importNetgenParams(hypo_file, myParams);
   fillHyp(myParams);
-
   MESSAGE("Meshing with netgen3d");
   int ret = Compute(myShape, *myMesh, myParams,
                       new_element_file,
index b31235bf95bf3dbc699280cd91ab7be7d9cef8a7..97037e779b55d9cf94b00008481b1ee2a527a067 100644 (file)
@@ -25,7 +25,9 @@
 //  Module : NETGEN
 //
 
+#include "NETGENPlugin_NETGEN_2D_SA.hxx"
 #include "NETGENPlugin_NETGEN_3D_SA.hxx"
+#include "NETGENPlugin_NETGEN_1D2D3D_SA.hxx"
 
 #include <stdio.h>
 #include <string.h>
@@ -78,18 +80,46 @@ int main(int argc, char *argv[]){
     element_orientation_file = "";
   if (new_element_file == "NONE")
     new_element_file = "";
-
+  int ret = 0;
   if (mesher=="NETGEN3D"){
     NETGENPlugin_NETGEN_3D_SA myplugin;
-    myplugin.run(input_mesh_file,
+    ret = myplugin.run(input_mesh_file,
              shape_file,
              hypo_file,
              element_orientation_file,
              new_element_file,
-             output_mesh_file);
-  } else {
+             output_mesh_file );
+  }
+  else if ( mesher=="NETGEN1D" || 
+            mesher=="NETGEN1D2D" || 
+            mesher=="NETGEN1D2D3D" ) 
+  {
+    NETGENPlugin_NETGEN_1D2D3D_SA myplugin;    
+    NETGENPlugin_Mesher::DIM DIM = mesher=="NETGEN1D" ? NETGENPlugin_Mesher::D1 
+                                    : ( mesher=="NETGEN1D2D" ? NETGENPlugin_Mesher::D2 
+                                    : NETGENPlugin_Mesher::D3 );
+
+    ret = myplugin.run(input_mesh_file,
+                  shape_file,
+                  hypo_file,
+                  element_orientation_file,
+                  new_element_file,
+                  output_mesh_file, 
+                  DIM );
+  } 
+  else if ( mesher=="NETGEN2D" ) 
+  {
+    NETGENPlugin_NETGEN_2D_SA myplugin;    
+    ret = myplugin.run(input_mesh_file,
+                  shape_file,
+                  hypo_file,
+                  element_orientation_file,
+                  new_element_file,
+                  output_mesh_file );
+  }  
+  else {
     std::cerr << "Unknown mesher:" << mesher << std::endl;
     return 1;
   }
-  return 0;
+  return ret;
 }
index 4e3d0cca889e942e2651f0e758d54decddf7a32e..e70b0346018756e2b8a6beef835c27b48b379642 100644 (file)
@@ -61,7 +61,7 @@ NETGENPlugin_SimpleHypothesis_2D::NETGENPlugin_SimpleHypothesis_2D (int
 //=============================================================================
 void NETGENPlugin_SimpleHypothesis_2D::SetNumberOfSegments(int nb)
 {
-  if ( nb < 1 )
+  if ( nb < 0 )
     throw SALOME_Exception("Number of segments must be positive");
   if (nb != _nbSegments)
   {
index 8dd6a5e9d7502e12da59dee50379c97ac3dc78fd..892bf3d42c4d14f7b13e18dbbdb1e7ffec3beebf 100644 (file)
@@ -37,6 +37,7 @@
 #include "NETGENPlugin_NETGEN_2D_i.hxx"
 #include "NETGENPlugin_NETGEN_3D_i.hxx"
 #include "NETGENPlugin_NETGEN_3D_Remote_i.hxx"
+#include "NETGENPlugin_NETGEN_2D_Remote_i.hxx"
 #include "NETGENPlugin_SimpleHypothesis_2D_i.hxx"
 #include "NETGENPlugin_SimpleHypothesis_3D_i.hxx"
 
@@ -66,6 +67,8 @@ extern "C"
       aCreator = new NETGENPlugin_Creator_i<NETGENPlugin_NETGEN_3D_i>;
     else if (strcmp(aHypName, "NETGEN_3D_Remote") == 0)
       aCreator = new NETGENPlugin_Creator_i<NETGENPlugin_NETGEN_3D_Remote_i>;
+    else if (strcmp(aHypName, "NETGEN_2D_Remote") == 0)
+      aCreator = new NETGENPlugin_Creator_i<NETGENPlugin_NETGEN_2D_Remote_i>;      
     else if (strcmp(aHypName, "NETGEN_2D") == 0)
       aCreator = new NETGENPlugin_Creator_i<NETGENPlugin_NETGEN_2D_i>;
     else if (strcmp(aHypName, "NETGEN_2D_ONLY") == 0)