src/NETGENPlugin/NETGENPlugin_NETGEN_2D_SA.cxx

   1 // Copyright (C) 2007-2024  CEA, EDF, OPEN CASCADE
   2 //
   3 // This library is free software; you can redistribute it and/or
   4 // modify it under the terms of the GNU Lesser General Public
   5 // License as published by the Free Software Foundation; either
   6 // version 2.1 of the License, or (at your option) any later version.
   7 //
   8 // This library is distributed in the hope that it will be useful,
   9 // but WITHOUT ANY WARRANTY; without even the implied warranty of
  10 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  11 // Lesser General Public License for more details.
  12 //
  13 // You should have received a copy of the GNU Lesser General Public
  14 // License along with this library; if not, write to the Free Software
  15 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
  16 //
  17 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
  18 //
  19
  20 // NETGENPlugin : C++ implementation
  21 // File   : NETGENPlugin_NETGEN_2D_SA.cxx
  22 // Author : Cesar Conopoima (OCC)
  23 // Date   : 23/10/2023
  24 // Project   : SALOME
  25 //=============================================================================
  26 //
  27 #include "NETGENPlugin_NETGEN_2D_SA.hxx"
  28 #include "NETGENPlugin_Hypothesis_2D.hxx"
  29 #include "NETGENPlugin_SimpleHypothesis_2D.hxx"
  30 #include "NETGENPlugin_Mesher.hxx"
  31
  32 #include <SMESHDS_Mesh.hxx>
  33 #include <SMESH_ControlsDef.hxx>
  34 #include <SMESH_Gen.hxx>
  35 #include <SMESH_Mesh.hxx>
  36 #include <SMESH_MesherHelper.hxx>
  37 #include <SMESH_DriverShape.hxx>
  38 #include <SMESH_DriverMesh.hxx>
  39 #include <StdMeshers_LengthFromEdges.hxx>
  40 #include <StdMeshers_MaxElementArea.hxx>
  41 #include <utilities.h>
  42
  43 #ifdef WIN32
  44 #include <filesystem>
  45 namespace fs = std::filesystem;
  46 #else
  47 #include <boost/filesystem.hpp>
  48 namespace fs = boost::filesystem;
  49 #endif
  50
  51 namespace nglib {
  52 #include <nglib.h>
  53 }
  54 #include <meshing.hpp>
  55
  56 using namespace nglib;
  57
  58 //=============================================================================
  59 /*!
  60  *
  61  */
  62 //=============================================================================
  63
  64 NETGENPlugin_NETGEN_2D_SA::NETGENPlugin_NETGEN_2D_SA()
  65   : NETGENPlugin_NETGEN_2D_ONLY(0, new SMESH_Gen())
  66 {
  67   _name = "NETGEN_2D_SA";
  68 }
  69
  70 //=============================================================================
  71 /*!
  72  *
  73  */
  74 //=============================================================================
  75
  76 NETGENPlugin_NETGEN_2D_SA::~NETGENPlugin_NETGEN_2D_SA()
  77 {
  78 }
  79
  80 /**
  81  * @brief Check presence and content of orientation file. Implemented for completness and future reference.
  82  *
  83  * @param element_orientation_file Binary file containing the orientation of surface elemnts
  84  * @return true, false
  85  */
  86 bool NETGENPlugin_NETGEN_2D_SA::checkOrientationFile( const std::string element_orientation_file )
  87 {
  88   if(element_orientation_file.empty()){
  89       MESSAGE("No element orientation file");
  90       return true;
  91     } else {
  92       MESSAGE("Reading from elements from file: " << element_orientation_file);
  93       // By construction the orientation file written by Remote version has a zero written to mark no need of orientation in 2D meshing
  94       int nbElement;
  95       std::ifstream df(element_orientation_file, ios::binary|ios::in);
  96       df.read((char*)&nbElement, sizeof(int));
  97       df.close();
  98       return (nbElement == 0);
  99     }
 100 }
 101
 102 /**
 103  * @brief fill plugin hypothesis from the netgen_params structure
 104  *
 105  * @param aParams the structure
 106  */
 107 void NETGENPlugin_NETGEN_2D_SA::fillHyp(const std::string param_file, netgen_params aParams)
 108 {
 109
 110   if( aParams.has_netgen_param && aParams.myType == hypoType::Hypo )
 111   {
 112     NETGENPlugin_Hypothesis_2D * hyp = new NETGENPlugin_Hypothesis_2D(0, GetGen());
 113     hyp->SetMaxSize(aParams.maxh);
 114     hyp->SetMinSize(aParams.minh);
 115     hyp->SetGrowthRate(aParams.grading);
 116     hyp->SetQuadAllowed(aParams.quad);
 117     hyp->SetFineness((NETGENPlugin_Hypothesis::Fineness)aParams.fineness);
 118     hyp->SetChordalErrorEnabled(aParams.chordalError);
 119     if(aParams.optimize){
 120       hyp->SetNbSurfOptSteps(aParams.optsteps2d);
 121     }
 122     _hypParameters = dynamic_cast< const NETGENPlugin_Hypothesis_2D *> (hyp);
 123   }
 124   else
 125   {
 126     NETGENPlugin_Hypothesis * hyp = new NETGENPlugin_Hypothesis(0, GetGen());
 127     if (param_file.find("lenghtfromedge") != std::string::npos )
 128     {
 129       _hypLengthFromEdges = dynamic_cast<const StdMeshers_LengthFromEdges*> (hyp);
 130     }
 131     else if ( param_file.find("maxarea") != std::string::npos )
 132     {
 133       StdMeshers_MaxElementArea * hyp = new StdMeshers_MaxElementArea(0, GetGen());
 134
 135       std::ifstream myfile(param_file);
 136       std::string line;
 137       double maxArea;
 138       std::getline(myfile, line);
 139       bool hashypothesis = std::stoi(line);
 140       if ( hashypothesis )
 141       {
 142         std::getline(myfile, line);
 143         maxArea = std::stod(line);
 144         hyp->SetMaxArea( maxArea );
 145         _hypMaxElementArea = static_cast<const StdMeshers_MaxElementArea*> (hyp);
 146       }
 147
 148       myfile.close();
 149     }
 150   }
 151 }
 152
 153
 154 /**
 155  * @brief Write a binary file containing information on the elements/nodes
 156  *        created by the mesher
 157  *
 158  * @param premeshedNodes map of the premeshed nodes of 1D elements
 159  * @param newNetgenCoordinates map of the coordinate of new netgen points created in the mesh
 160  * @param newNetgenElements map of the element triangulation
 161  * @return true if there are some error
 162  */
 163 bool NETGENPlugin_NETGEN_2D_SA::fillNewElementFile( std::string new_element_file,
 164                                                     const int numberOfGlobalPremeshedNodes,
 165                                                     std::map<int,const SMDS_MeshNode*>& premeshedNodes,
 166                                                     std::map<int,std::vector<double>>& newNetgenCoordinates,
 167                                                     std::map<int,std::vector<smIdType>>& newNetgenElements )
 168 {
 169   MESSAGE("Writting new elements")
 170
 171   int NetgenNbOfNodes     = premeshedNodes.size();
 172   int NetgenNbOfNodesNew  = numberOfGlobalPremeshedNodes + newNetgenCoordinates.size();
 173   int NetgenNbOfTriangles = newNetgenElements.size();
 174   bool isOK = ( NetgenNbOfTriangles > 0 );
 175   if ( isOK && !new_element_file.empty() )
 176   {
 177     int    NetgenElement[3];
 178     std::ofstream df(new_element_file, ios::out|ios::binary);
 179     // Writing nodevec (correspondence netgen numbering mesh numbering)
 180     // Number of nodes
 181     df.write((char*) &numberOfGlobalPremeshedNodes, sizeof(int));
 182     df.write((char*) &NetgenNbOfNodes, sizeof(int));
 183     df.write((char*) &NetgenNbOfNodesNew, sizeof(int));
 184     for (auto k : premeshedNodes )
 185       df.write((char*) &k.first, sizeof(int));
 186
 187     // Writing info on new points
 188     for (auto k : newNetgenCoordinates )
 189       df.write((char*) k.second.data(), sizeof(double)*3 );
 190
 191     // create triangles
 192     df.write((char*) &NetgenNbOfTriangles, sizeof(int));
 193     for ( int elemIndex = 1; elemIndex <= NetgenNbOfTriangles; ++elemIndex )
 194     {
 195       int nodes = newNetgenElements[ elemIndex ].size();
 196       for (int i = 0; i < nodes; i++)
 197         NetgenElement[ i ] = (int) newNetgenElements[ elemIndex ][ i ];
 198
 199       df.write((char*) &NetgenElement, sizeof(int)* nodes );
 200     }
 201
 202     df.close();
 203   }
 204   return false;
 205 }
 206
 207 /**
 208  * @brief Compute the mesh based on the
 209  *
 210  * @param input_mesh_file Mesh file (containing 2D elements)
 211  * @param shape_file Shape file (BREP or STEP format)
 212  * @param hypo_file Ascii file containing the netgen parameters
 213  * @param element_orientation_file Binary file containing the orientation of surface elemnts
 214  * @param new_element_file output file containing info the elements created by the mesher
 215  * @param output_mesh_file output mesh file (if empty it will not be created)
 216  * @return true, false
 217  */
 218 bool NETGENPlugin_NETGEN_2D_SA::Compute( SMESH_Mesh& aMesh, TopoDS_Shape &aShape, std::string new_element_file )
 219 {
 220   // Nodes on edge are double because each face is treat one by one, so seam edges are iterated twice
 221   vector< const SMDS_MeshNode* > nodeVec;
 222   NETGENPlugin_Mesher aMesher( &aMesh, aShape, false /*isVolume=*/ );
 223   NETGENPlugin_NetgenLibWrapper ngLib;
 224   ngLib._isComputeOk = false;
 225
 226   fs::path netgen_log_file = fs::path(new_element_file).remove_filename() / fs::path("NETGEN.out");
 227   MESSAGE("netgen ouput"<<netgen_log_file.string());
 228   ngLib.setOutputFile(netgen_log_file.string());
 229
 230   netgen::OCCGeometry occgeoComm;
 231   NETGENPlugin_NETGEN_2D_ONLY::SetParameteres( aMesh, aShape, aMesher, ngLib._ngMesh, occgeoComm, false /*submesh is not supported*/ );
 232
 233   std::map<int,const SMDS_MeshNode*> premeshedNodes;
 234   std::map<int,std::vector<double>> newNetgenCoordinates;
 235   std::map<int,std::vector<smIdType>> newNetgenElements;
 236   const int numberOfTotalPremeshedNodes = aMesh.NbNodes();
 237   bool compute = NETGENPlugin_NETGEN_2D_ONLY::MapSegmentsToEdges( aMesh, aShape, ngLib, nodeVec,
 238                                                                   premeshedNodes, newNetgenCoordinates,
 239                                                                   newNetgenElements );
 240
 241   compute = fillNewElementFile(new_element_file,
 242                                 numberOfTotalPremeshedNodes,
 243                                 premeshedNodes,
 244                                 newNetgenCoordinates,
 245                                 newNetgenElements);
 246   return compute;
 247 }
 248
 249 /**
 250  * @brief Running the mesher on the given files
 251  *
 252  * @param input_mesh_file Mesh file (containing 2D elements)
 253  * @param shape_file Shape file (BREP or STEP format)
 254  * @param hypo_file Ascii file containing the netgen parameters
 255  * @param element_orientation_file Binary file containing the orientation of surface elemnts
 256  * @param new_element_file output file containing info the elements created by the mesher
 257  * @param output_mesh_file output mesh file (if empty it will not be created)
 258  * @return int
 259  */
 260 int NETGENPlugin_NETGEN_2D_SA::run(const std::string input_mesh_file,
 261                                     const std::string shape_file,
 262                                     const std::string hypo_file,
 263                                     const std::string element_orientation_file,
 264                                     const std::string new_element_file,
 265                                     const std::string output_mesh_file)
 266 {
 267
 268   std::unique_ptr<SMESH_Mesh> myMesh(_gen->CreateMesh(false));
 269
 270   SMESH_DriverMesh::importMesh(input_mesh_file, *myMesh);
 271   // Importing shape
 272   TopoDS_Shape myShape;
 273   SMESH_DriverShape::importShape(shape_file, myShape);
 274   // Importing hypothesis
 275   netgen_params myParams;
 276   importNetgenParams(hypo_file, myParams);
 277   fillHyp(hypo_file,myParams);
 278   MESSAGE("Meshing with netgen2d");
 279   int ret = 1;
 280   if ( checkOrientationFile(element_orientation_file) )
 281   {
 282     ret = (int) Compute( *myMesh, myShape, new_element_file );
 283
 284     if(ret){
 285       std::cerr << "Meshing failed" << std::endl;
 286       return ret;
 287     }
 288
 289     if(!output_mesh_file.empty()){
 290       std::string meshName = "MESH";
 291       SMESH_DriverMesh::exportMesh(output_mesh_file, *myMesh, meshName);
 292     }
 293   }
 294   else
 295     std::cerr << "For NETGENPlugin_NETGEN_2D_SA orientation file should be market with 0 or be empty!" << std::endl;
 296
 297   return ret;
 298 }