ForkMaster/netgenplugin
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

Creating Remote and SA NETGEN_3D plugins

Browse Source
This commit is contained in:
Yoann Audouin 2022-09-21 14:21:30 +02:00
parent a52ccf7016
commit 5e49461800
15 changed files with 1034 additions and 1380 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

7
idl/NETGENPlugin_Algorithm.idl
View File

@ -44,6 +44,13 @@ module NETGENPlugin
{
};
/*!
* NETGENPlugin_NETGEN_3D: interface of "Remote Tetrahedron (Netgen)" algorithm
*/
interface NETGENPlugin_NETGEN_3D_Remote : NETGENPlugin::NETGENPlugin_NETGEN_3D
{
};
/*!
* NETGENPlugin_NETGEN_2D: interface of "Netgen 1D-2D" algorithm
*/

8
src/NETGENPlugin/CMakeLists.txt
View File

@ -99,7 +99,9 @@ SET(NETGENEngine_HEADERS
NETGENPlugin_Remesher_2D.hxx
NETGENPlugin_Defs.hxx
NETGENPlugin_DriverParam.hxx
NETGENPlugin_Runner.hxx
NETGENPlugin_NETGEN_3D_SA.hxx
NETGENPlugin_NETGEN_3D_Remote.hxx
NETGENPlugin_NETGEN_3D_Remote_i.hxx
)
# --- sources ---
@ -127,8 +129,10 @@ SET(NETGENEngine_SOURCES
NETGENPlugin_SimpleHypothesis_3D_i.cxx
NETGENPlugin_Remesher_2D.cxx
NETGENPlugin_i.cxx
NETGENPlugin_Runner.cxx
NETGENPlugin_DriverParam.cxx
NETGENPlugin_NETGEN_3D_SA.cxx
NETGENPlugin_NETGEN_3D_Remote.cxx
NETGENPlugin_NETGEN_3D_Remote_i.cxx
)
SET(NetgenRunner_SOURCES

32
src/NETGENPlugin/NETGENPluginBuilder.py
View File

@ -42,6 +42,7 @@ NETGEN_VERSION_MAJOR = NETGENPlugin.NETGEN_VERSION_MAJOR
## Algorithm type: Netgen tetrahedron 3D algorithm, see NETGEN_3D_Algorithm
NETGEN_3D = "NETGEN_3D"
NETGEN_3D_Remote = "NETGEN_3D_Remote"
## Algorithm type: Netgen tetrahedron 1D-2D-3D algorithm, see NETGEN_1D2D3D_Algorithm
NETGEN_1D2D3D = "NETGEN_2D3D"
## Algorithm type: Netgen triangle 1D-2D algorithm, see NETGEN_1D2D_Algorithm
@ -470,6 +471,37 @@ class NETGEN_3D_Algorithm(NETGEN_Algorithm):
pass # end of NETGEN_3D_Algorithm class
## Tetrahedron 3D algorithm
#
# It can be created by calling smeshBuilder.Mesh.Tetrahedron() or smeshBuilder.Mesh.Tetrahedron( smeshBuilder.NETGEN, geom=0 )
#
# This algorithm generates only 3D (volumes) elements for given geometrical shape
# and, in contrast to NETGEN_1D2D3D_Algorithm class, should be used in conjunction
# with other 1D and 2D meshing algorithms.
class NETGEN_3D_Remote_Algorithm(NETGEN_3D_Algorithm):
## type of algorithm used with helper function in smeshBuilder.Mesh class
# @internal
algoType = NETGEN_3D_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 tetrahedron 3D algorithm for 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):
NETGEN_3D_Algorithm.__init__(self, mesh, geom)
pass
pass # end of NETGEN_3D_Remote_Algorithm class
## Triangle (helper) 1D-2D algorithm
#

280
src/NETGENPlugin/NETGENPlugin_NETGEN_3D.cxx
View File

@ -32,6 +32,7 @@
#include "NETGENPlugin_Hypothesis.hxx"
// TODO: remove use of netgen_param
#include "NETGENPlugin_DriverParam.hxx"
#include <SMDS_MeshElement.hxx>
@ -48,8 +49,6 @@
#include <StdMeshers_QuadToTriaAdaptor.hxx>
#include <StdMeshers_ViscousLayers.hxx>
#include <SMESH_subMesh.hxx>
#include <SMESH_DriverShape.hxx>
#include <SMESH_DriverMesh.hxx>
#include <BRepGProp.hxx>
@ -69,12 +68,7 @@
#include <vector>
#include <map>
#include <QString>
#include <QProcess>
#include <cstdlib>
#include <boost/filesystem.hpp>
namespace fs = boost::filesystem;
/*
Netgen include files
@ -200,273 +194,6 @@ bool NETGENPlugin_NETGEN_3D::CheckHypothesis (SMESH_Mesh& aMesh,
}
void NETGENPlugin_NETGEN_3D::fillParameters(const NETGENPlugin_Hypothesis* hyp, netgen_params &aParams)
{
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();
#else
// const char*
aParams.meshsizefilename = hyp->GetMeshSizeFile();
#endif
#ifdef NETGEN_V6
aParams.closeedgefac = 2;
#else
aParams.closeedgefac = 0;
#endif
}
// write in a binary file the orientation for each 2D element of the mesh
void NETGENPlugin_NETGEN_3D::exportElementOrientation(SMESH_Mesh& aMesh,
const TopoDS_Shape& aShape,
netgen_params& aParams,
const std::string output_file)
{
SMESH_MesherHelper helper(aMesh);
NETGENPlugin_Internals internals( aMesh, aShape, /*is3D=*/true );
SMESH_ProxyMesh::Ptr proxyMesh( new SMESH_ProxyMesh( aMesh ));
std::map<vtkIdType, bool> elemOrientation;
for ( TopExp_Explorer exFa( aShape, TopAbs_FACE ); exFa.More(); exFa.Next())
{
const TopoDS_Shape& aShapeFace = exFa.Current();
int faceID = aMesh.GetMeshDS()->ShapeToIndex( aShapeFace );
bool isInternalFace = internals.isInternalShape( faceID );
bool isRev = false;
if ( !isInternalFace &&
helper.NbAncestors(aShapeFace, aMesh, aShape.ShapeType()) > 1 )
// IsReversedSubMesh() can work wrong on strongly curved faces,
// so we use it as less as possible
isRev = helper.IsReversedSubMesh( TopoDS::Face( aShapeFace ));
const SMESHDS_SubMesh * aSubMeshDSFace = proxyMesh->GetSubMesh( aShapeFace );
if ( !aSubMeshDSFace ) continue;
SMDS_ElemIteratorPtr iteratorElem = aSubMeshDSFace->GetElements();
if ( aParams._quadraticMesh &&
dynamic_cast< const SMESH_ProxyMesh::SubMesh*>( aSubMeshDSFace ))
{
// add medium nodes of proxy triangles to helper (#16843)
while ( iteratorElem->more() )
helper.AddTLinks( static_cast< const SMDS_MeshFace* >( iteratorElem->next() ));
iteratorElem = aSubMeshDSFace->GetElements();
}
while ( iteratorElem->more() ) // loop on elements on a geom face
{
// check mesh face
const SMDS_MeshElement* elem = iteratorElem->next();
if ( !elem )
error( COMPERR_BAD_INPUT_MESH, "Null element encounters");
if ( elem->NbCornerNodes() != 3 )
error( COMPERR_BAD_INPUT_MESH, "Not triangle element encounters");
elemOrientation[elem->GetID()] = isRev;
} // loop on elements on a face
} // loop on faces of a SOLID or SHELL
{
std::ofstream df(output_file, ios::out|ios::binary);
int size=elemOrientation.size();
df.write((char*)&size, sizeof(int));
for(auto const& [id, orient]:elemOrientation){
df.write((char*)&id, sizeof(vtkIdType));
df.write((char*)&orient, sizeof(bool));
}
}
}
int NETGENPlugin_NETGEN_3D::RemoteCompute(SMESH_Mesh& aMesh,
const TopoDS_Shape& aShape)
{
aMesh.Lock();
auto time0 = std::chrono::high_resolution_clock::now();
SMESH_Hypothesis::Hypothesis_Status hypStatus;
CheckHypothesis(aMesh, aShape, hypStatus);
auto time1 = std::chrono::high_resolution_clock::now();
auto elapsed = std::chrono::duration_cast<std::chrono::nanoseconds>(time1-time0);
std::cout << "Time for check_hypo: " << elapsed.count() * 1e-9 << std::endl;
// Temporary folder for run
fs::path tmp_folder = aMesh.tmp_folder / fs::unique_path(fs::path("Volume-%%%%-%%%%"));
fs::create_directories(tmp_folder);
// Using MESH2D generated after all triangles where created.
fs::path mesh_file=aMesh.tmp_folder / fs::path("Mesh2D.med");
fs::path element_orientation_file=tmp_folder / fs::path("element_orientation.dat");
fs::path new_element_file=tmp_folder / fs::path("new_elements.dat");
fs::path tmp_mesh_file=tmp_folder / fs::path("tmp_mesh.med");
// 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("netgen3d_param.txt");
fs::path log_file=tmp_folder / fs::path("run.log");
fs::path cmd_file=tmp_folder / fs::path("cmd.log");
//TODO: Handle variable mesh_name
std::string mesh_name = "Maillage_1";
//Writing Shape
exportShape(shape_file.string(), aShape);
auto time2 = std::chrono::high_resolution_clock::now();
elapsed = std::chrono::duration_cast<std::chrono::nanoseconds>(time2-time1);
std::cout << "Time for exportShape: " << elapsed.count() * 1e-9 << std::endl;
//Writing hypo
netgen_params aParams;
fillParameters(_hypParameters, aParams);
exportNetgenParams(param_file.string(), aParams);
auto time3 = std::chrono::high_resolution_clock::now();
elapsed = std::chrono::duration_cast<std::chrono::nanoseconds>(time3-time2);
std::cout << "Time for fill+export param: " << elapsed.count() * 1e-9 << std::endl;
// Exporting element orientation
exportElementOrientation(aMesh, aShape, aParams, element_orientation_file.string());
auto time4 = std::chrono::high_resolution_clock::now();
elapsed = std::chrono::duration_cast<std::chrono::nanoseconds>(time4-time3);
std::cout << "Time for exportElemOrient: " << elapsed.count() * 1e-9 << std::endl;
aMesh.Unlock();
// Calling run_mesher
// TODO: check if we need to handle the .exe for windows
std::string cmd;
fs::path run_mesher_exe =
fs::path(std::getenv("NETGENPLUGIN_ROOT_DIR"))/
fs::path("bin")/
fs::path("salome")/
fs::path("NETGENPlugin_Runner");
cmd = run_mesher_exe.string() +
" NETGEN3D " + mesh_file.string() + " "
+ shape_file.string() + " "
+ param_file.string() + " "
+ element_orientation_file.string() + " "
+ std::to_string(aMesh.GetMesherNbThreads()) + " "
+ new_element_file.string() + " "
+ "NONE";
// Writing command in log
{
std::ofstream flog(cmd_file.string());
flog << cmd << endl;
flog << endl;
}
//std::cout << "Running command: " << std::endl;
//std::cout << cmd << std::endl;
// Building arguments for QProcess
QString program = run_mesher_exe.c_str();
QStringList arguments;
arguments << "NETGEN3D";
arguments << mesh_file.c_str();
arguments << shape_file.c_str();
arguments << param_file.c_str();
arguments << element_orientation_file.c_str();
arguments << std::to_string(aMesh.GetMesherNbThreads()).c_str();
arguments << new_element_file.c_str();
arguments << "NONE";
QString out_file = log_file.c_str();
QProcess myProcess;
myProcess.setStandardOutputFile(out_file);
myProcess.start(program, arguments);
myProcess.waitForFinished();
int ret = myProcess.exitStatus();
auto time5 = std::chrono::high_resolution_clock::now();
elapsed = std::chrono::duration_cast<std::chrono::nanoseconds>(time5-time4);
std::cout << "Time for exec of run_mesher: " << elapsed.count() * 1e-9 << std::endl;
if(ret != 0){
// Run crahed
std::cout << "Issue with command: " << std::endl;
std::cout << "See log for more details: " << log_file.string() << std::endl;
std::cout << cmd << std::endl;
return false;
}
aMesh.Lock();
{
std::ifstream df(new_element_file.string(), ios::binary);
int Netgen_NbOfNodes;
int Netgen_NbOfNodesNew;
int Netgen_NbOfTetra;
double Netgen_point[3];
int Netgen_tetrahedron[4];
int nodeID;
SMESH_MesherHelper helper(aMesh);
// This function is mandatory for setElementsOnShape to work
helper.IsQuadraticSubMesh(aShape);
helper.SetElementsOnShape( true );
// Number of nodes in intial mesh
df.read((char*) &Netgen_NbOfNodes, sizeof(int));
// Number of nodes added by netgen
df.read((char*) &Netgen_NbOfNodesNew, sizeof(int));
// Filling nodevec (correspondence netgen numbering mesh numbering)
vector< const SMDS_MeshNode* > nodeVec ( Netgen_NbOfNodesNew + 1 );
//vector<int> nodeTmpVec ( Netgen_NbOfNodesNew + 1 );
SMESHDS_Mesh * meshDS = helper.GetMeshDS();
for (int nodeIndex = 1 ; nodeIndex <= Netgen_NbOfNodes; ++nodeIndex )
{
//Id of the point
df.read((char*) &nodeID, sizeof(int));
nodeVec.at(nodeIndex) = meshDS->FindNode(nodeID);
}
// Add new points and update nodeVec
for (int nodeIndex = Netgen_NbOfNodes +1 ; nodeIndex <= Netgen_NbOfNodesNew; ++nodeIndex )
{
df.read((char *) &Netgen_point, sizeof(double)*3);
nodeVec.at(nodeIndex) = helper.AddNode(Netgen_point[0],
Netgen_point[1],
Netgen_point[2]);
}
// Add tetrahedrons
df.read((char*) &Netgen_NbOfTetra, sizeof(int));
for ( int elemIndex = 1; elemIndex <= Netgen_NbOfTetra; ++elemIndex )
{
df.read((char*) &Netgen_tetrahedron, sizeof(int)*4);
helper.AddVolume(
nodeVec.at( Netgen_tetrahedron[0] ),
nodeVec.at( Netgen_tetrahedron[1] ),
nodeVec.at( Netgen_tetrahedron[2] ),
nodeVec.at( Netgen_tetrahedron[3] ));
}
}
auto time7 = std::chrono::high_resolution_clock::now();
elapsed = std::chrono::duration_cast<std::chrono::nanoseconds>(time7-time5);
std::cout << "Time for exec of add_in_mesh: " << elapsed.count() * 1e-9 << std::endl;
aMesh.Unlock();
return true;
}
//=============================================================================
/*!
@ -479,7 +206,7 @@ int NETGENPlugin_NETGEN_3D::RemoteCompute(SMESH_Mesh& aMesh,
*
*/
bool getSurfaceElements(
bool NETGENPlugin_NETGEN_3D::getSurfaceElements(
SMESH_Mesh& aMesh,
const TopoDS_Shape& aShape,
SMESH_ProxyMesh::Ptr proxyMesh,
@ -851,9 +578,6 @@ bool NETGENPlugin_NETGEN_3D::Compute(
SMESH_Mesh& aMesh,
const TopoDS_Shape& aShape)
{
if(aMesh.IsParallel())
return RemoteCompute(aMesh, aShape);
// vector of nodes in which node index == netgen ID
vector< const SMDS_MeshNode* > nodeVec;
NETGENPlugin_NetgenLibWrapper ngLib;

27
src/NETGENPlugin/NETGENPlugin_NETGEN_3D.hxx
View File

@ -38,6 +38,7 @@
#include "Utils_SALOME_Exception.hxx"
#include <vector>
#include <tuple>
class StdMeshers_ViscousLayers;
class StdMeshers_MaxElementVolume;
@ -70,7 +71,7 @@ class NETGENPLUGIN_EXPORT NETGENPlugin_NETGEN_3D: public SMESH_3D_Algo
const TopoDS_Shape& aShape,
MapShapeNbElems& aResMap);
static bool computeFillNgMesh(
bool computeFillNgMesh(
SMESH_Mesh& aMesh,
const TopoDS_Shape& aShape,
std::vector< const SMDS_MeshNode* > &nodeVec,
@ -79,7 +80,7 @@ class NETGENPLUGIN_EXPORT NETGENPlugin_NETGEN_3D: public SMESH_3D_Algo
netgen_params &aParams,
int &Netgen_NbOfNodes);
static bool computePrepareParam(
bool computePrepareParam(
SMESH_Mesh& aMesh,
NETGENPlugin_NetgenLibWrapper &ngLib,
netgen::OCCGeometry &occgeo,
@ -87,7 +88,7 @@ class NETGENPLUGIN_EXPORT NETGENPlugin_NETGEN_3D: public SMESH_3D_Algo
netgen_params &aParams,
int &endWith);
static bool computeRunMesher(
bool computeRunMesher(
netgen::OCCGeometry &occgeo,
std::vector< const SMDS_MeshNode* > &nodeVec,
netgen::Mesh* ngMesh,
@ -95,7 +96,7 @@ class NETGENPLUGIN_EXPORT NETGENPlugin_NETGEN_3D: public SMESH_3D_Algo
netgen_params &aParams,
int &startWith, int &endWith);
static bool computeFillMesh(
bool computeFillMesh(
std::vector< const SMDS_MeshNode* > &nodeVec,
NETGENPlugin_NetgenLibWrapper &ngLib,
SMESH_MesherHelper &helper,
@ -103,18 +104,14 @@ class NETGENPLUGIN_EXPORT NETGENPlugin_NETGEN_3D: public SMESH_3D_Algo
protected:
void exportElementOrientation(SMESH_Mesh& aMesh,
virtual bool getSurfaceElements(
SMESH_Mesh& aMesh,
const TopoDS_Shape& aShape,
netgen_params& aParams,
const std::string output_file);
void fillParameters(const NETGENPlugin_Hypothesis* hyp,
netgen_params &aParams);
int RemoteCompute(SMESH_Mesh& aMesh,
const TopoDS_Shape& aShape);
SMESH_ProxyMesh::Ptr proxyMesh,
NETGENPlugin_Internals &internals,
SMESH_MesherHelper &helper,
netgen_params &aParams,
std::map<const SMDS_MeshElement*, std::tuple<bool, bool>>& listElements);
bool compute(SMESH_Mesh& mesh,
SMESH_MesherHelper& helper,

351
src/NETGENPlugin/NETGENPlugin_NETGEN_3D_Remote.cxx Normal file
View File

@ -0,0 +1,351 @@
// Copyright (C) 2007-2021 CEA/DEN, EDF R&D, 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_3D_Remote.cxx
// Created : lundi 19 Septembre 2022
// Author : Yoann AUDOUIN (CEA)
// Project : SALOME
//=============================================================================
//
//
#include "NETGENPlugin_NETGEN_3D_Remote.hxx"
#include "NETGENPlugin_NETGEN_3D.hxx"
#include "NETGENPlugin_DriverParam.hxx"
#include "NETGENPlugin_Hypothesis.hxx"
#include <SMESH_Gen.hxx>
#include <SMESH_Mesh.hxx>
#include <SMESH_MesherHelper.hxx>
#include <SMESH_DriverShape.hxx>
#include <SMESH_DriverMesh.hxx>
#include <SMESHDS_Mesh.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;
//=============================================================================
/*!
*
*/
//=============================================================================
NETGENPlugin_NETGEN_3D_Remote::NETGENPlugin_NETGEN_3D_Remote(int hypId, SMESH_Gen * gen)
: NETGENPlugin_NETGEN_3D(hypId, gen)
{
_name = "NETGEN_3D_Remote";
}
//=============================================================================
/*!
*
*/
//=============================================================================
NETGENPlugin_NETGEN_3D_Remote::~NETGENPlugin_NETGEN_3D_Remote()
{
}
void NETGENPlugin_NETGEN_3D_Remote::fillParameters(const NETGENPlugin_Hypothesis* hyp, netgen_params &aParams)
{
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();
#else
// const char*
aParams.meshsizefilename = hyp->GetMeshSizeFile();
#endif
#ifdef NETGEN_V6
aParams.closeedgefac = 2;
#else
aParams.closeedgefac = 0;
#endif
}
// write in a binary file the orientation for each 2D element of the mesh
void NETGENPlugin_NETGEN_3D_Remote::exportElementOrientation(SMESH_Mesh& aMesh,
const TopoDS_Shape& aShape,
netgen_params& aParams,
const std::string output_file)
{
SMESH_MesherHelper helper(aMesh);
NETGENPlugin_Internals internals( aMesh, aShape, /*is3D=*/true );
SMESH_ProxyMesh::Ptr proxyMesh( new SMESH_ProxyMesh( aMesh ));
std::map<vtkIdType, bool> elemOrientation;
for ( TopExp_Explorer exFa( aShape, TopAbs_FACE ); exFa.More(); exFa.Next())
{
const TopoDS_Shape& aShapeFace = exFa.Current();
int faceID = aMesh.GetMeshDS()->ShapeToIndex( aShapeFace );
bool isInternalFace = internals.isInternalShape( faceID );
bool isRev = false;
if ( !isInternalFace &&
helper.NbAncestors(aShapeFace, aMesh, aShape.ShapeType()) > 1 )
// IsReversedSubMesh() can work wrong on strongly curved faces,
// so we use it as less as possible
isRev = helper.IsReversedSubMesh( TopoDS::Face( aShapeFace ));
const SMESHDS_SubMesh * aSubMeshDSFace = proxyMesh->GetSubMesh( aShapeFace );
if ( !aSubMeshDSFace ) continue;
SMDS_ElemIteratorPtr iteratorElem = aSubMeshDSFace->GetElements();
if ( aParams._quadraticMesh &&
dynamic_cast< const SMESH_ProxyMesh::SubMesh*>( aSubMeshDSFace ))
{
// add medium nodes of proxy triangles to helper (#16843)
while ( iteratorElem->more() )
helper.AddTLinks( static_cast< const SMDS_MeshFace* >( iteratorElem->next() ));
iteratorElem = aSubMeshDSFace->GetElements();
}
while ( iteratorElem->more() ) // loop on elements on a geom face
{
// check mesh face
const SMDS_MeshElement* elem = iteratorElem->next();
if ( !elem )
error( COMPERR_BAD_INPUT_MESH, "Null element encounters");
if ( elem->NbCornerNodes() != 3 )
error( COMPERR_BAD_INPUT_MESH, "Not triangle element encounters");
elemOrientation[elem->GetID()] = isRev;
} // loop on elements on a face
} // loop on faces of a SOLID or SHELL
{
std::ofstream df(output_file, ios::out|ios::binary);
int size=elemOrientation.size();
df.write((char*)&size, sizeof(int));
for(auto const& [id, orient]:elemOrientation){
df.write((char*)&id, sizeof(vtkIdType));
df.write((char*)&orient, sizeof(bool));
}
}
}
bool NETGENPlugin_NETGEN_3D_Remote::Compute(SMESH_Mesh& aMesh,
const TopoDS_Shape& aShape)
{
aMesh.Lock();
SMESH_Hypothesis::Hypothesis_Status hypStatus;
NETGENPlugin_NETGEN_3D::CheckHypothesis(aMesh, aShape, hypStatus);
// Temporary folder for run
fs::path tmp_folder = aMesh.tmp_folder / fs::unique_path(fs::path("Volume-%%%%-%%%%"));
fs::create_directories(tmp_folder);
// Using MESH2D generated after all triangles where created.
fs::path mesh_file=aMesh.tmp_folder / fs::path("Mesh2D.med");
fs::path element_orientation_file=tmp_folder / fs::path("element_orientation.dat");
fs::path new_element_file=tmp_folder / fs::path("new_elements.dat");
fs::path tmp_mesh_file=tmp_folder / fs::path("tmp_mesh.med");
// 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("netgen3d_param.txt");
fs::path log_file=tmp_folder / fs::path("run.log");
fs::path cmd_file=tmp_folder / fs::path("cmd.log");
//TODO: Handle variable mesh_name
std::string mesh_name = "Maillage_1";
//Writing Shape
exportShape(shape_file.string(), aShape);
//Writing hypo
netgen_params aParams;
std::cout << _hypParameters << std::endl;
fillParameters(_hypParameters, aParams);
exportNetgenParams(param_file.string(), aParams);
// Exporting element orientation
exportElementOrientation(aMesh, aShape, aParams, element_orientation_file.string());
aMesh.Unlock();
// Calling run_mesher
// TODO: check if we need to handle the .exe for windows
std::string cmd;
fs::path run_mesher_exe =
fs::path(std::getenv("NETGENPLUGIN_ROOT_DIR"))/
fs::path("bin")/
fs::path("salome")/
fs::path("NETGENPlugin_Runner");
cmd = run_mesher_exe.string() +
" NETGEN3D " + mesh_file.string() + " "
+ shape_file.string() + " "
+ param_file.string() + " "
+ element_orientation_file.string() + " "
+ std::to_string(aMesh.GetMesherNbThreads()) + " "
+ new_element_file.string() + " "
+ "NONE";
// Writing command in log
{
std::ofstream flog(cmd_file.string());
flog << cmd << endl;
flog << endl;
}
//std::cout << "Running command: " << std::endl;
//std::cout << cmd << std::endl;
// Building arguments for QProcess
QString program = run_mesher_exe.c_str();
QStringList arguments;
arguments << "NETGEN3D";
arguments << mesh_file.c_str();
arguments << shape_file.c_str();
arguments << param_file.c_str();
arguments << element_orientation_file.c_str();
arguments << std::to_string(aMesh.GetMesherNbThreads()).c_str();
arguments << new_element_file.c_str();
arguments << "NONE";
QString out_file = log_file.c_str();
QProcess myProcess;
myProcess.setStandardOutputFile(out_file);
myProcess.start(program, arguments);
myProcess.waitForFinished();
int ret = myProcess.exitStatus();
if(ret != 0){
// Run crahed
std::cout << "Issue with command: " << std::endl;
std::cout << "See log for more details: " << log_file.string() << std::endl;
std::cout << cmd << std::endl;
return false;
}
aMesh.Lock();
{
std::ifstream df(new_element_file.string(), ios::binary);
int Netgen_NbOfNodes;
int Netgen_NbOfNodesNew;
int Netgen_NbOfTetra;
double Netgen_point[3];
int Netgen_tetrahedron[4];
int nodeID;
SMESH_MesherHelper helper(aMesh);
// This function is mandatory for setElementsOnShape to work
helper.IsQuadraticSubMesh(aShape);
helper.SetElementsOnShape( true );
// Number of nodes in intial mesh
df.read((char*) &Netgen_NbOfNodes, sizeof(int));
// Number of nodes added by netgen
df.read((char*) &Netgen_NbOfNodesNew, sizeof(int));
// Filling nodevec (correspondence netgen numbering mesh numbering)
vector< const SMDS_MeshNode* > nodeVec ( Netgen_NbOfNodesNew + 1 );
//vector<int> nodeTmpVec ( Netgen_NbOfNodesNew + 1 );
SMESHDS_Mesh * meshDS = helper.GetMeshDS();
for (int nodeIndex = 1 ; nodeIndex <= Netgen_NbOfNodes; ++nodeIndex )
{
//Id of the point
df.read((char*) &nodeID, sizeof(int));
nodeVec.at(nodeIndex) = meshDS->FindNode(nodeID);
}
// Add new points and update nodeVec
for (int nodeIndex = Netgen_NbOfNodes +1 ; nodeIndex <= Netgen_NbOfNodesNew; ++nodeIndex )
{
df.read((char *) &Netgen_point, sizeof(double)*3);
nodeVec.at(nodeIndex) = helper.AddNode(Netgen_point[0],
Netgen_point[1],
Netgen_point[2]);
}
// Add tetrahedrons
df.read((char*) &Netgen_NbOfTetra, sizeof(int));
for ( int elemIndex = 1; elemIndex <= Netgen_NbOfTetra; ++elemIndex )
{
df.read((char*) &Netgen_tetrahedron, sizeof(int)*4);
helper.AddVolume(
nodeVec.at( Netgen_tetrahedron[0] ),
nodeVec.at( Netgen_tetrahedron[1] ),
nodeVec.at( Netgen_tetrahedron[2] ),
nodeVec.at( Netgen_tetrahedron[3] ));
}
}
aMesh.Unlock();
return true;
}

74
src/NETGENPlugin/NETGENPlugin_NETGEN_3D_Remote.hxx Normal file
View File

@ -0,0 +1,74 @@
// Copyright (C) 2007-2021 CEA/DEN, EDF R&D, 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_3D_Remote.hxx
// Created : lundi 19 Septembre 2022
// Author : Yoann AUDOUIN (EDF)
// Project : SALOME
//=============================================================================
//
#ifndef _NETGENPlugin_NETGEN_3D_REMOTE_HXX_
#define _NETGENPlugin_NETGEN_3D_REMOTE_HXX_
#include "NETGENPlugin_NETGEN_3D.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_3D_Remote: public NETGENPlugin_NETGEN_3D
{
public:
NETGENPlugin_NETGEN_3D_Remote(int hypId, SMESH_Gen* gen);
virtual ~NETGENPlugin_NETGEN_3D_Remote();
// Function whould not be used with remote Computing
bool CheckHypothesis (SMESH_Mesh& aMesh,
const TopoDS_Shape& aShape,
Hypothesis_Status& aStatus){aStatus = HYP_OK;return true;};
bool Compute(SMESH_Mesh& aMesh,
const TopoDS_Shape& aShape);
protected:
void exportElementOrientation(SMESH_Mesh& aMesh,
const TopoDS_Shape& aShape,
netgen_params& aParams,
const std::string output_file);
void fillParameters(const NETGENPlugin_Hypothesis* hyp,
netgen_params &aParams);
};
#endif

80
src/NETGENPlugin/NETGENPlugin_NETGEN_3D_Remote_i.cxx Normal file
View File

@ -0,0 +1,80 @@
// Copyright (C) 2007-2021 CEA/DEN, EDF R&D, 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_3D_Remote_i.cxx
// Author : Yoann AUDOUIN (EDF)
// Module : NETGENPlugin
// $Header$
//
#include "NETGENPlugin_NETGEN_3D_Remote_i.hxx"
#include "SMESH_Gen.hxx"
#include "Utils_CorbaException.hxx"
#include "utilities.h"
using namespace std;
//=============================================================================
/*!
* NETGENPlugin_NETGEN_3D_Remote_i::NETGENPlugin_NETGEN_3D_Remote_i
*
* Constructor
*/
//=============================================================================
NETGENPlugin_NETGEN_3D_Remote_i::NETGENPlugin_NETGEN_3D_Remote_i( PortableServer::POA_ptr thePOA,
::SMESH_Gen* theGenImpl )
: SALOME::GenericObj_i( thePOA ),
SMESH_Hypothesis_i( thePOA ),
SMESH_Algo_i( thePOA ),
SMESH_3D_Algo_i( thePOA )
{
myBaseImpl = new ::NETGENPlugin_NETGEN_3D_Remote( theGenImpl->GetANewId(),
theGenImpl );
}
//=============================================================================
/*!
* NETGENPlugin_NETGEN_3D_Remote_i::~NETGENPlugin_NETGEN_3D_Remote_i
*
* Destructor
*/
//=============================================================================
NETGENPlugin_NETGEN_3D_Remote_i::~NETGENPlugin_NETGEN_3D_Remote_i()
{
}
//=============================================================================
/*!
* NETGENPlugin_NETGEN_3D_Remote_i::GetImpl
*
* Get implementation
*/
//=============================================================================
::NETGENPlugin_NETGEN_3D_Remote* NETGENPlugin_NETGEN_3D_Remote_i::GetImpl()
{
return ( ::NETGENPlugin_NETGEN_3D_Remote* )myBaseImpl;
}

58
src/NETGENPlugin/NETGENPlugin_NETGEN_3D_Remote_i.hxx Normal file
View File

@ -0,0 +1,58 @@
// Copyright (C) 2007-2021 CEA/DEN, EDF R&D, 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_3D_Remote_i.hxx
// Author : Yoann AUDOUIN (EDF)
// Module : NETGENPlugin
// $Header$
//
#ifndef _NETGENPlugin_NETGEN_3D_REMOTE_I_HXX_
#define _NETGENPlugin_NETGEN_3D_REMOTE_I_HXX_
#include "NETGENPlugin_Defs.hxx"
#include <SALOMEconfig.h>
#include CORBA_SERVER_HEADER(NETGENPlugin_Algorithm)
#include "SMESH_3D_Algo_i.hxx"
#include "NETGENPlugin_NETGEN_3D_Remote.hxx"
// ======================================================
// NETGEN 3d algorithm
// ======================================================
class NETGENPLUGIN_EXPORT NETGENPlugin_NETGEN_3D_Remote_i:
public virtual POA_NETGENPlugin::NETGENPlugin_NETGEN_3D_Remote,
public virtual SMESH_3D_Algo_i
{
public:
// Constructor
NETGENPlugin_NETGEN_3D_Remote_i( PortableServer::POA_ptr thePOA,
::SMESH_Gen* theGenImpl );
// Destructor
virtual ~NETGENPlugin_NETGEN_3D_Remote_i();
// Get implementation
::NETGENPlugin_NETGEN_3D_Remote* GetImpl();
};
#endif

314
src/NETGENPlugin/NETGENPlugin_NETGEN_3D_SA.cxx Normal file
View File

@ -0,0 +1,314 @@
// Copyright (C) 2007-2021 CEA/DEN, EDF R&D, 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_3D_SA.cxx
// Created : lundi 19 Septembre 2022
// Author : Yoann AUDOUIN (CEA)
// Project : SALOME
//=============================================================================
//
//
#include "NETGENPlugin_NETGEN_3D_SA.hxx"
#include "NETGENPlugin_DriverParam.hxx"
#include "NETGENPlugin_Hypothesis.hxx"
#include <SMESH_Gen.hxx>
#include <SMESH_Mesh.hxx>
#include <SMESH_MesherHelper.hxx>
#include <SMESH_DriverShape.hxx>
#include <SMESH_DriverMesh.hxx>
#include <SMESHDS_Mesh.hxx>
#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;
//=============================================================================
/*!
*
*/
//=============================================================================
NETGENPlugin_NETGEN_3D_SA::NETGENPlugin_NETGEN_3D_SA()
: NETGENPlugin_NETGEN_3D(0, new SMESH_Gen())
{
_name = "NETGEN_3D_SA";
}
//=============================================================================
/*!
*
*/
//=============================================================================
NETGENPlugin_NETGEN_3D_SA::~NETGENPlugin_NETGEN_3D_SA()
{
}
bool NETGENPlugin_NETGEN_3D_SA::computeFillNewElementFile(
std::vector< const SMDS_MeshNode* > &nodeVec,
NETGENPlugin_NetgenLibWrapper &ngLib,
std::string new_element_file,
int &Netgen_NbOfNodes
)
{
Ng_Mesh* Netgen_mesh = ngLib.ngMesh();
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() )
{
std::ofstream df(new_element_file, ios::out|ios::binary);
double Netgen_point[3];
int Netgen_tetrahedron[4];
// Writing nodevec (correspondence netgen numbering mesh numbering)
// Number of nodes
df.write((char*) &Netgen_NbOfNodes, sizeof(int));
df.write((char*) &Netgen_NbOfNodesNew, sizeof(int));
for (int nodeIndex = 1 ; nodeIndex <= Netgen_NbOfNodes; ++nodeIndex )
{
//Id of the point
int id = nodeVec.at(nodeIndex)->GetID();
df.write((char*) &id, sizeof(int));
}
// Writing info on new points
for (int nodeIndex = Netgen_NbOfNodes +1 ; nodeIndex <= Netgen_NbOfNodesNew; ++nodeIndex )
{
Ng_GetPoint(Netgen_mesh, nodeIndex, Netgen_point );
// Coordinates of the point
df.write((char *) &Netgen_point, sizeof(double)*3);
}
// create tetrahedrons
df.write((char*) &Netgen_NbOfTetra, sizeof(int));
for ( int elemIndex = 1; elemIndex <= Netgen_NbOfTetra; ++elemIndex )
{
Ng_GetVolumeElement(Netgen_mesh, elemIndex, Netgen_tetrahedron);
df.write((char*) &Netgen_tetrahedron, sizeof(int)*4);
}
}
return false;
}
bool NETGENPlugin_NETGEN_3D_SA::Compute(TopoDS_Shape &aShape, SMESH_Mesh& aMesh, netgen_params& aParams,
std::string new_element_file, std::string element_orientation_file,
bool output_mesh)
{
// vector of nodes in which node index == netgen ID
vector< const SMDS_MeshNode* > nodeVec;
NETGENPlugin_NetgenLibWrapper ngLib;
SMESH_MesherHelper helper(aMesh);
int startWith = netgen::MESHCONST_MESHVOLUME;
int endWith = netgen::MESHCONST_OPTVOLUME;
int Netgen_NbOfNodes=0;
bool ret;
ret = NETGENPlugin_NETGEN_3D::computeFillNgMesh(aMesh, aShape, nodeVec, ngLib, helper, aParams, Netgen_NbOfNodes);
if(ret)
return error( aParams._error, aParams._comment);
netgen::OCCGeometry occgeo;
NETGENPlugin_NETGEN_3D::computePrepareParam(aMesh, ngLib, occgeo, helper, aParams, endWith);
ret = NETGENPlugin_NETGEN_3D::computeRunMesher(occgeo, nodeVec, ngLib._ngMesh, ngLib, aParams, startWith, endWith);
if(ret){
if(aParams._error)
return error(aParams._error, aParams._comment);
error(aParams._comment);
return true;
}
computeFillNewElementFile(nodeVec, ngLib, new_element_file, Netgen_NbOfNodes);
if(output_mesh)
NETGENPlugin_NETGEN_3D::computeFillMesh(nodeVec, ngLib, helper, Netgen_NbOfNodes);
return false;
}
int NETGENPlugin_NETGEN_3D_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,
int nbThreads)
{
_element_orientation_file = element_orientation_file;
// Importing mesh
SMESH_Gen gen;
std::unique_ptr<SMESH_Mesh> myMesh(gen.CreateMesh(false));
//TODO: To define
std::string mesh_name = "Maillage_1";
importMesh(input_mesh_file, *myMesh, mesh_name);
// Importing shape
TopoDS_Shape myShape;
importShape(shape_file, myShape);
// Importing hypothesis
netgen_params myParams;
importNetgenParams(hypo_file, myParams, &gen);
// Setting number of threads for netgen
myParams.nbThreads = nbThreads;
std::cout << "Meshing with netgen3d" << std::endl;
int ret = Compute(myShape, *myMesh, myParams,
new_element_file, element_orientation_file,
!output_mesh_file.empty());
if(ret){
std::cout << "Meshing failed" << std::endl;
return ret;
}
if(!output_mesh_file.empty()){
exportMesh(output_mesh_file, *myMesh, mesh_name);
}
return ret;
}
bool NETGENPlugin_NETGEN_3D_SA::getSurfaceElements(
SMESH_Mesh& aMesh,
const TopoDS_Shape& aShape,
SMESH_ProxyMesh::Ptr proxyMesh,
NETGENPlugin_Internals &internals,
SMESH_MesherHelper &helper,
netgen_params &aParams,
std::map<const SMDS_MeshElement*, tuple<bool, bool>>& listElements
)
{
SMESHDS_Mesh* meshDS = aMesh.GetMeshDS();
// Get list of elements + their orientation from element_orientation file
std::map<vtkIdType, bool> elemOrientation;
{
// Setting all element orientation to false if there no element orientation file
if(_element_orientation_file.empty()){
std::cout << "No element orientation file" << std::endl;
SMDS_ElemIteratorPtr iteratorElem = meshDS->elementsIterator(SMDSAbs_Face);
while ( iteratorElem->more() ) // loop on elements on a geom face
{
// check mesh face
const SMDS_MeshElement* elem = iteratorElem->next();
elemOrientation[elem->GetID()] = false;
}
} else {
std::cout << "Reading from elements from file: " << _element_orientation_file << std::endl;
std::ifstream df(_element_orientation_file, ios::binary|ios::in);
int nbElement;
bool orient;
// Warning of the use of vtkIdType (I had issue when run_mesher was compiled with internal vtk) and salome not
// Sizeof was the same but how he othered the type was different
// Maybe using another type (uint64_t) instead would be better
vtkIdType id;
df.read((char*)&nbElement, sizeof(int));
for(int ielem=0;ielem<nbElement;++ielem){
df.read((char*) &id, sizeof(vtkIdType));
df.read((char*) &orient, sizeof(bool));
elemOrientation[id] = orient;
}
}
}
// Adding elements from Mesh
SMDS_ElemIteratorPtr iteratorElem = meshDS->elementsIterator(SMDSAbs_Face);
bool isRev;
bool isInternalFace = false;
bool isIn;
while ( iteratorElem->more() ) // loop on elements on a geom face
{
// check mesh face
const SMDS_MeshElement* elem = iteratorElem->next();
if ( !elem ){
aParams._error = COMPERR_BAD_INPUT_MESH;
aParams._comment = "Null element encounters";
return true;
}
if ( elem->NbCornerNodes() != 3 ){
aParams._error = COMPERR_BAD_INPUT_MESH;
aParams._comment = "Not triangle element encounters";
return true;
}
// Keeping only element that are in the element orientation file
isIn = elemOrientation.count(elem->GetID())==1;
if(!isIn)
continue;
// Get orientation
// Netgen requires that all the triangle point outside
isRev = elemOrientation[elem->GetID()];
listElements[elem] = tuple(isRev, false);
}
return false;
}

89
src/NETGENPlugin/NETGENPlugin_NETGEN_3D_SA.hxx Normal file
View File

@ -0,0 +1,89 @@
// Copyright (C) 2007-2021 CEA/DEN, EDF R&D, 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_3D_SA.hxx
// Created : lundi 19 Septembre 2022
// Author : Yoann AUDOUIN (EDF)
// Project : SALOME
//=============================================================================
//
#ifndef _NETGENPlugin_NETGEN_3D_SA_HXX_
#define _NETGENPlugin_NETGEN_3D_SA_HXX_
#include "NETGENPlugin_NETGEN_3D.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_3D_SA: public NETGENPlugin_NETGEN_3D
{
public:
NETGENPlugin_NETGEN_3D_SA();
virtual ~NETGENPlugin_NETGEN_3D_SA();
bool Compute(TopoDS_Shape &aShape, SMESH_Mesh& aMesh, netgen_params& aParams,
std::string new_element_file, std::string element_orientation_file,
bool output_mesh);
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,
int nbThreads);
protected:
bool computeFillNewElementFile(
std::vector< const SMDS_MeshNode* > &nodeVec,
NETGENPlugin_NetgenLibWrapper &ngLib,
std::string new_element_file,
int &Netgen_NbOfNodes);
bool getSurfaceElements(
SMESH_Mesh& aMesh,
const TopoDS_Shape& aShape,
SMESH_ProxyMesh::Ptr proxyMesh,
NETGENPlugin_Internals &internals,
SMESH_MesherHelper &helper,
netgen_params &aParams,
std::map<const SMDS_MeshElement*, tuple<bool, bool>>& listElements
);
std::string _element_orientation_file="";
};
#endif

955
src/NETGENPlugin/NETGENPlugin_Runner.cxx
View File

@ -1,955 +0,0 @@
// Copyright (C) 2007-2021 CEA/DEN, EDF R&D, 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_Runner.cxx
// Author : Yoann AUDOUIN, EDF
// Module : SMESH
//
#include "NETGENPlugin_Runner.hxx"
#include "NETGENPlugin_NETGEN_3D.hxx"
#include "NETGENPlugin_DriverParam.hxx"
#include <fstream>
#include <vector>
#include <boost/filesystem.hpp>
namespace fs = boost::filesystem;
#include <chrono>
// SMESH include
#include <SMESH_Mesh.hxx>
#include <SMESH_subMesh.hxx>
#include <SMESH_Gen.hxx>
#include <SMESH_Algo.hxx>
#include <SMESHDS_Mesh.hxx>
#include <SMESH_ControlsDef.hxx>
#include <SMESH_Comment.hxx>
#include <SMESH_ComputeError.hxx>
#include <SMESH_MesherHelper.hxx>
#include <StdMeshers_MaxElementVolume.hxx>
#include <StdMeshers_QuadToTriaAdaptor.hxx>
#include <StdMeshers_ViscousLayers.hxx>
#include <StdMeshers_ViscousLayers2D.hxx>
#include <SMESH_DriverShape.hxx>
#include <SMESH_DriverMesh.hxx>
// NETGENPlugin
// #include <NETGENPlugin_Mesher.hxx>
// #include <NETGENPlugin_Hypothesis.hxx>
#include "NETGENPlugin_Mesher.hxx"
#include "NETGENPlugin_Hypothesis.hxx"
// OCC include
#include <TopoDS.hxx>
#include <BRepClass3d_SolidClassifier.hxx>
#include <GProp_GProps.hxx>
#include <BRepGProp.hxx>
#include <Standard_Failure.hxx>
#include <Standard_ErrorHandler.hxx>
/*
Netgen include files
*/
#ifndef OCCGEOMETRY
#define OCCGEOMETRY
#endif
#include <occgeom.hpp>
#include <meshing.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;
NETGENPLUGIN_DLL_HEADER
extern bool merge_solids;
#ifdef NETGEN_V5
extern void OCCSetLocalMeshSize(OCCGeometry & geom, Mesh & mesh);
#endif
}
using namespace nglib;
int error(int error_type, std::string msg)
{
std::cerr << msg << std::endl;
return error_type;
};
int error(const SMESH_Comment& comment)
{
return error(1, "SMESH_Comment error: "+comment);
};
/**
* @brief Set the netgen parameters
*
* @param aParams Internal structure of parameters
* @param mparams Netgen strcuture of parameters
*/
void set_netgen_parameters(netgen_params& aParams)
{
// Default parameters
#ifdef NETGEN_V6
//netgen::mparam.nthreads = std::thread::hardware_concurrency();
netgen::mparam.nthreads = aParams.nbThreads;
netgen::mparam.parallel_meshing = aParams.nbThreads > 1;
if ( getenv( "SALOME_NETGEN_DISABLE_MULTITHREADING" ))
{
netgen::mparam.nthreads = 1;
netgen::mparam.parallel_meshing = false;
}
#endif
// Initialize global NETGEN parameters:
netgen::mparam.maxh = aParams.maxh;
netgen::mparam.minh = aParams.minh;
netgen::mparam.segmentsperedge = aParams.segmentsperedge;
netgen::mparam.grading = aParams.grading;
netgen::mparam.curvaturesafety = aParams.curvaturesafety;
netgen::mparam.secondorder = aParams.secondorder;
netgen::mparam.quad = aParams.quad;
netgen::mparam.uselocalh = aParams.uselocalh;
netgen::merge_solids = aParams.merge_solids;
netgen::mparam.optsteps2d = aParams.optsteps2d;
netgen::mparam.optsteps3d = aParams.optsteps3d;
netgen::mparam.elsizeweight = aParams.elsizeweight;
netgen::mparam.opterrpow = aParams.opterrpow;
netgen::mparam.delaunay = aParams.delaunay;
netgen::mparam.checkoverlap = aParams.checkoverlap;
netgen::mparam.checkchartboundary = aParams.checkchartboundary;
#ifdef NETGEN_V6
// std::string
netgen::mparam.meshsizefilename = aParams.meshsizefilename;
netgen::mparam.closeedgefac = aParams.closeedgefac;
#else
// const char*
netgen::mparam.meshsizefilename= aParams.meshsizefilename.empty() ? 0 : aParams.meshsizefilename.c_str();
#endif
}
/**
* @brief compute mesh with netgen3d
*
* @param input_mesh_file Input Mesh file
* @param shape_file Shape file
* @param hypo_file Parameter file
* @param new_element_file Binary file containing new nodes and new element info
* @param output_mesh If true will export mesh into output_mesh_file
* @param output_mesh_file Output Mesh file
*
* @return error code
*/
int netgen3d(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,
int nbThreads)
{
// Importing mesh
SMESH_Gen gen;
std::unique_ptr<SMESH_Mesh> myMesh(gen.CreateMesh(false));
//TODO: To define
std::string mesh_name = "Maillage_1";
importMesh(input_mesh_file, *myMesh, mesh_name);
// Importing shape
TopoDS_Shape myShape;
importShape(shape_file, myShape);
// Importing hypothesis
netgen_params myParams;
importNetgenParams(hypo_file, myParams, &gen);
// Setting number of threads for netgen
myParams.nbThreads = nbThreads;
std::cout << "Meshing with netgen3d" << std::endl;
int ret = netgen3dInternal(myShape, *myMesh, myParams,
new_element_file, element_orientation_file,
!output_mesh_file.empty());
if(ret){
std::cout << "Meshing failed" << std::endl;
return ret;
}
if(!output_mesh_file.empty()){
exportMesh(output_mesh_file, *myMesh, mesh_name);
}
return ret;
}
bool getSurfaceElements(
SMESH_Mesh& aMesh,
const TopoDS_Shape& aShape,
SMESH_ProxyMesh::Ptr proxyMesh,
NETGENPlugin_Internals &internals,
SMESH_MesherHelper &helper,
netgen_params &aParams,
std::string element_orientation_file,
std::map<const SMDS_MeshElement*, tuple<bool, bool>>& listElements
)
{
SMESHDS_Mesh* meshDS = aMesh.GetMeshDS();
// Get list of elements + their orientation from element_orientation file
std::map<vtkIdType, bool> elemOrientation;
{
// Setting all element orientation to false if there no element orientation file
if(element_orientation_file.empty()){
SMDS_ElemIteratorPtr iteratorElem = meshDS->elementsIterator(SMDSAbs_Face);
while ( iteratorElem->more() ) // loop on elements on a geom face
{
// check mesh face
const SMDS_MeshElement* elem = iteratorElem->next();
elemOrientation[elem->GetID()] = false;
}
} else {
std::ifstream df(element_orientation_file, ios::binary|ios::in);
int nbElement;
bool orient;
// Warning of the use of vtkIdType (I had issue when run_mesher was compiled with internal vtk) and salome not
// Sizeof was the same but how he othered the type was different
// Maybe using another type (uint64_t) instead would be better
vtkIdType id;
df.read((char*)&nbElement, sizeof(int));
for(int ielem=0;ielem<nbElement;++ielem){
df.read((char*) &id, sizeof(vtkIdType));
df.read((char*) &orient, sizeof(bool));
elemOrientation[id] = orient;
}
}
}
// Adding elements from Mesh
SMDS_ElemIteratorPtr iteratorElem = meshDS->elementsIterator(SMDSAbs_Face);
bool isRev;
bool isInternalFace = false;
bool isIn;
while ( iteratorElem->more() ) // loop on elements on a geom face
{
// check mesh face
const SMDS_MeshElement* elem = iteratorElem->next();
if ( !elem ){
aParams._error = COMPERR_BAD_INPUT_MESH;
aParams._comment = "Null element encounters";
return true;
}
if ( elem->NbCornerNodes() != 3 ){
aParams._error = COMPERR_BAD_INPUT_MESH;
aParams._comment = "Not triangle element encounters";
return true;
}
// Keeping only element that are in the element orientation file
isIn = elemOrientation.count(elem->GetID())==1;
if(!isIn)
continue;
// Get orientation
// Netgen requires that all the triangle point outside
isRev = elemOrientation[elem->GetID()];
listElements[elem] = tuple(isRev, false);
}
return false;
}
bool mycomputeFillNgMesh(
SMESH_Mesh& aMesh,
const TopoDS_Shape& aShape,
std::vector< const SMDS_MeshNode* > &nodeVec,
NETGENPlugin_NetgenLibWrapper &ngLib,
SMESH_MesherHelper &helper,
netgen_params &aParams,
std::string element_orientation_file,
int &Netgen_NbOfNodes)
{
netgen::multithread.terminate = 0;
netgen::multithread.task = "Volume meshing";
aParams._progressByTic = -1.;
SMESHDS_Mesh* meshDS = aMesh.GetMeshDS();
aParams._quadraticMesh = helper.IsQuadraticSubMesh(aShape);
helper.SetElementsOnShape( true );
Netgen_NbOfNodes = 0;
double Netgen_point[3];
int Netgen_triangle[3];
Ng_Mesh * Netgen_mesh = (Ng_Mesh*)ngLib._ngMesh;
{
const int invalid_ID = -1;
SMESH::Controls::Area areaControl;
SMESH::Controls::TSequenceOfXYZ nodesCoords;
// maps nodes to ng ID
typedef map< const SMDS_MeshNode*, int, TIDCompare > TNodeToIDMap;
typedef TNodeToIDMap::value_type TN2ID;
TNodeToIDMap nodeToNetgenID;
// find internal shapes
NETGENPlugin_Internals internals( aMesh, aShape, /*is3D=*/true );
// ---------------------------------
// Feed the Netgen with surface mesh
// ---------------------------------
bool isRev=false;
bool isInternalFace=false;
SMESH_ProxyMesh::Ptr proxyMesh( new SMESH_ProxyMesh( aMesh ));
if ( aParams._viscousLayersHyp )
{
netgen::multithread.percent = 3;
proxyMesh = aParams._viscousLayersHyp->Compute( aMesh, aShape );
if ( !proxyMesh )
return false;
}
if ( aMesh.NbQuadrangles() > 0 )
{
netgen::multithread.percent = 6;
StdMeshers_QuadToTriaAdaptor* Adaptor = new StdMeshers_QuadToTriaAdaptor;
Adaptor->Compute(aMesh,aShape,proxyMesh.get());
proxyMesh.reset( Adaptor );
}
std::map<const SMDS_MeshElement*, tuple<bool, bool>> listElements;
bool ret = getSurfaceElements(aMesh, aShape, proxyMesh, internals, helper, aParams, element_orientation_file, listElements);
if(ret)
return ret;
for ( auto const& [elem, info] : listElements ) // loop on elements on a geom face
{
isRev = get<0>(info);
isInternalFace = get<1>(info);
// Add nodes of triangles and triangles them-selves to netgen mesh
// add three nodes of triangle
bool hasDegen = false;
for ( int iN = 0; iN < 3; ++iN )
{
const SMDS_MeshNode* node = elem->GetNode( iN );
const int shapeID = node->getshapeId();
if ( node->GetPosition()->GetTypeOfPosition() == SMDS_TOP_EDGE &&
helper.IsDegenShape( shapeID ))
{
// ignore all nodes on degeneraged edge and use node on its vertex instead
TopoDS_Shape vertex = TopoDS_Iterator( meshDS->IndexToShape( shapeID )).Value();
node = SMESH_Algo::VertexNode( TopoDS::Vertex( vertex ), meshDS );
hasDegen = true;
}
int& ngID = nodeToNetgenID.insert(TN2ID( node, invalid_ID )).first->second;
if ( ngID == invalid_ID )
{
ngID = ++Netgen_NbOfNodes;
Netgen_point [ 0 ] = node->X();
Netgen_point [ 1 ] = node->Y();
Netgen_point [ 2 ] = node->Z();
Ng_AddPoint(Netgen_mesh, Netgen_point);
}
Netgen_triangle[ isRev ? 2-iN : iN ] = ngID;
}
// add triangle
if ( hasDegen && (Netgen_triangle[0] == Netgen_triangle[1] ||
Netgen_triangle[0] == Netgen_triangle[2] ||
Netgen_triangle[2] == Netgen_triangle[1] ))
continue;
Ng_AddSurfaceElement(Netgen_mesh, NG_TRIG, Netgen_triangle);
if ( isInternalFace && !proxyMesh->IsTemporary( elem ))
{
swap( Netgen_triangle[1], Netgen_triangle[2] );
Ng_AddSurfaceElement(Netgen_mesh, NG_TRIG, Netgen_triangle);
}
} // loop on elements on a face
// insert old nodes into nodeVec
nodeVec.resize( nodeToNetgenID.size() + 1, 0 );
TNodeToIDMap::iterator n_id = nodeToNetgenID.begin();
for ( ; n_id != nodeToNetgenID.end(); ++n_id )
nodeVec[ n_id->second ] = n_id->first;
nodeToNetgenID.clear();
if ( internals.hasInternalVertexInSolid() )
{
netgen::OCCGeometry occgeo;
NETGENPlugin_Mesher::AddIntVerticesInSolids( occgeo,
(netgen::Mesh&) *Netgen_mesh,
nodeVec,
internals);
}
}
Netgen_NbOfNodes = Ng_GetNP( Netgen_mesh );
return false;
}
bool mycomputeFillNewElementFile(
std::vector< const SMDS_MeshNode* > &nodeVec,
NETGENPlugin_NetgenLibWrapper &ngLib,
std::string new_element_file,
int &Netgen_NbOfNodes
)
{
Ng_Mesh* Netgen_mesh = ngLib.ngMesh();
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() )
{
std::ofstream df(new_element_file, ios::out|ios::binary);
double Netgen_point[3];
int Netgen_tetrahedron[4];
// Writing nodevec (correspondence netgen numbering mesh numbering)
// Number of nodes
df.write((char*) &Netgen_NbOfNodes, sizeof(int));
df.write((char*) &Netgen_NbOfNodesNew, sizeof(int));
for (int nodeIndex = 1 ; nodeIndex <= Netgen_NbOfNodes; ++nodeIndex )
{
//Id of the point
int id = nodeVec.at(nodeIndex)->GetID();
df.write((char*) &id, sizeof(int));
}
// Writing info on new points
for (int nodeIndex = Netgen_NbOfNodes +1 ; nodeIndex <= Netgen_NbOfNodesNew; ++nodeIndex )
{
Ng_GetPoint(Netgen_mesh, nodeIndex, Netgen_point );
// Coordinates of the point
df.write((char *) &Netgen_point, sizeof(double)*3);
}
// create tetrahedrons
df.write((char*) &Netgen_NbOfTetra, sizeof(int));
for ( int elemIndex = 1; elemIndex <= Netgen_NbOfTetra; ++elemIndex )
{
Ng_GetVolumeElement(Netgen_mesh, elemIndex, Netgen_tetrahedron);
df.write((char*) &Netgen_tetrahedron, sizeof(int)*4);
}
}
return false;
}
/**
* @brief Compute aShape within aMesh using netgen3d
*
* @param aShape the shape
* @param aMesh the mesh
* @param aParams the netgen parameters
* @param new_element_file file containing data on the new point/tetra added by netgen
* @param element_orientation_file file containing data on the orientation of each element to add to netgen
* @param output_mesh if true add element created by netgen into aMesh
*
* @return error code
*/
int netgen3dInternal(TopoDS_Shape &aShape, SMESH_Mesh& aMesh, netgen_params& aParams,
std::string new_element_file, std::string element_orientation_file,
bool output_mesh)
{
// vector of nodes in which node index == netgen ID
vector< const SMDS_MeshNode* > nodeVec;
NETGENPlugin_NetgenLibWrapper ngLib;
SMESH_MesherHelper helper(aMesh);
int startWith = netgen::MESHCONST_MESHVOLUME;
int endWith = netgen::MESHCONST_OPTVOLUME;
int Netgen_NbOfNodes=0;
bool ret;
ret = mycomputeFillNgMesh(aMesh, aShape, nodeVec, ngLib, helper, aParams, element_orientation_file, Netgen_NbOfNodes);
if(ret)
return error( aParams._error, aParams._comment);
netgen::OCCGeometry occgeo;
NETGENPlugin_NETGEN_3D::computePrepareParam(aMesh, ngLib, occgeo, helper, aParams, endWith);
ret = NETGENPlugin_NETGEN_3D::computeRunMesher(occgeo, nodeVec, ngLib._ngMesh, ngLib, aParams, startWith, endWith);
if(ret){
if(aParams._error)
return error(aParams._error, aParams._comment);
error(aParams._comment);
return true;
}
mycomputeFillNewElementFile(nodeVec, ngLib, new_element_file, Netgen_NbOfNodes);
if(output_mesh)
NETGENPlugin_NETGEN_3D::computeFillMesh(nodeVec, ngLib, helper, Netgen_NbOfNodes);
return false;
}
/**
* @brief compute mesh with netgen2d
*
* @param input_mesh_file Input Mesh file
* @param shape_file Shape file
* @param hypo_file Parameter file
* @param new_element_file Binary file containing new nodes and new element info
* @param output_mesh If true will export mesh into output_mesh_file
* @param output_mesh_file Output Mesh file
*
* @return error code
*/
int netgen2d(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)
{
// Importing mesh
SMESH_Gen gen;
std::unique_ptr<SMESH_Mesh> myMesh(gen.CreateMesh(false));
//TODO: To define
std::string mesh_name = "Maillage_1";
importMesh(input_mesh_file, *myMesh, mesh_name);
// Importing shape
TopoDS_Shape myShape;
importShape(shape_file, myShape);
// Importing hypothesis
netgen_params myParams;
importNetgenParams(hypo_file, myParams, &gen);
std::cout << "Meshing with netgen3d" << std::endl;
int ret = netgen2dInternal(myShape, *myMesh, myParams,
new_element_file, element_orientation_file,
!output_mesh_file.empty());
if(!ret){
std::cout << "Meshing failed" << std::endl;
return ret;
}
if(!output_mesh_file.empty())
exportMesh(output_mesh_file, *myMesh, mesh_name);
return ret;
}
// TODO: Not working properly
/**
* @brief Compute aShape within aMesh using netgen2d
*
* @param aShape the shape
* @param aMesh the mesh
* @param aParams the netgen parameters
* @param new_element_file file containing data on the new point/tetra added by netgen
*
* @return error code
*/
int netgen2dInternal(TopoDS_Shape &aShape, SMESH_Mesh& aMesh, netgen_params& aParams,
std::string new_element_file, std::string element_orientation_file,
bool output_mesh)
{
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;
std::map<vtkIdType, bool> elemOrientation;
typedef map< const SMDS_MeshNode*, int, TIDCompare > TNodeToIDMap;
typedef TNodeToIDMap::value_type TN2ID;
const int invalid_ID = -1;
int Netgen_NbOfNodes=0;
double Netgen_point[3];
int Netgen_segment[2];
int Netgen_triangle[3];
// min / max sizes are set as follows:
// if ( _hypParameters )
// min and max are defined by the user
// else if ( aParams.has_LengthFromEdges_hyp )
// 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);
set_netgen_parameters( aParams );
const bool toOptimize = aParams.optimize;
if ( aParams.has_maxelementvolume_hyp )
{
netgen::mparam.maxh = sqrt( 2. * aParams.maxElementVolume / sqrt(3.0) );
}
netgen::mparam.quad = aParams.quad;
// local size is common for all FACEs in aShape?
const bool isCommonLocalSize = ( !aParams.has_LengthFromEdges_hyp && !aParams.has_maxelementvolume_hyp && netgen::mparam.uselocalh );
const bool isDefaultHyp = ( !aParams.has_LengthFromEdges_hyp && !aParams.has_maxelementvolume_hyp && !aParams.has_netgen_param );
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
// minh, face_maxh, grading and curvaturesafety; find minh if not set by the user
if ( !aParams.has_netgen_param || netgen::mparam.minh < DBL_MIN )
{
if ( !aParams.has_netgen_param )
netgen::mparam.maxh = occgeoComm.GetBoundingBox().Diam() / 3.;
netgen::mparam.minh = aMesher.GetDefaultMinSize( aShape, netgen::mparam.maxh );
}
// set local size depending on curvature and NOT closeness of EDGEs
#ifdef NETGEN_V6
const double factor = 2; //netgen::occparam.resthcloseedgefac;
#else
const double factor = netgen::occparam.resthcloseedgefac;
netgen::occparam.resthcloseedgeenable = false;
netgen::occparam.resthcloseedgefac = 1.0 + netgen::mparam.grading;
#endif
occgeoComm.face_maxh = netgen::mparam.maxh;
#ifdef NETGEN_V6
netgen::OCCParameters occparam;
netgen::OCCSetLocalMeshSize( occgeoComm, *ngMeshes[0], netgen::mparam, occparam );
#else
netgen::OCCSetLocalMeshSize( occgeoComm, *ngMeshes[0] );
#endif
occgeoComm.emap.Clear();
occgeoComm.vmap.Clear();
// Reading list of element to integrate into netgen mesh
{
std::ifstream df(element_orientation_file, ios::in|ios::binary);
int nbElement;
vtkIdType id;
bool orient;
df.read((char*)&nbElement, sizeof(int));
for(int ielem=0;ielem<nbElement;++ielem){
df.read((char*) &id, sizeof(vtkIdType));
df.read((char*) &orient, sizeof(bool));
elemOrientation[id] = orient;
}
}
bool isIn;
// set local size according to size of existing segments
SMDS_ElemIteratorPtr iteratorElem = meshDS->elementsIterator(SMDSAbs_Edge);
while ( iteratorElem->more() ) // loop on elements on a geom face
{
const SMDS_MeshElement* seg = iteratorElem->next();
// Keeping only element that are in the element orientation file
isIn = elemOrientation.count(seg->GetID())==1;
if(!isIn)
continue;
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[0]->RestrictLocalH( pi, factor * ( n1 - n2 ).Modulus() );
}
// set local size defined on shapes
aMesher.SetLocalSize( occgeoComm, *ngMeshes[0] );
aMesher.SetLocalSizeForChordalError( occgeoComm, *ngMeshes[0] );
try {
ngMeshes[0]->LoadLocalMeshSize( netgen::mparam.meshsizefilename );
} catch (netgen::NgException & ex) {
return error( COMPERR_BAD_PARMETERS, ex.What() );
}
}
netgen::mparam.uselocalh = toOptimize; // restore as it is used at surface optimization
// ==================
// Loop on all FACEs
// ==================
vector< const SMDS_MeshNode* > nodeVec;
// prepare occgeom
netgen::OCCGeometry occgeom;
occgeom.shape = aShape;
occgeom.fmap.Add( aShape );
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;
// -------------------------
// Fill netgen mesh
// -------------------------
// maps nodes to ng ID
// MESHCONST_ANALYSE step may lead to a failure, so we make an attempt
// w/o MESHCONST_ANALYSE at the second loop
int err = 0;
enum { LOC_SIZE, NO_LOC_SIZE };
int iLoop = isCommonLocalSize ? 0 : 1;
int faceID = occgeom.fmap.FindIndex(aShape);
int solidID = 0;
for ( ; iLoop < 2; iLoop++ )
{
//bool isMESHCONST_ANALYSE = false;
//TODO: check how to replace that
//InitComputeError();
netgen::Mesh * ngMesh = ngMeshes[ iLoop ];
ngMesh->DeleteMesh();
if ( iLoop == NO_LOC_SIZE )
{
ngMesh->SetGlobalH ( netgen::mparam.maxh );
ngMesh->SetMinimalH( netgen::mparam.minh );
netgen::Box<3> bb = occgeom.GetBoundingBox();
bb.Increase (bb.Diam()/10);
ngMesh->SetLocalH (bb.PMin(), bb.PMax(), netgen::mparam.grading);
aMesher.SetLocalSize( occgeom, *ngMesh );
aMesher.SetLocalSizeForChordalError( occgeoComm, *ngMesh );
try {
ngMesh->LoadLocalMeshSize( netgen::mparam.meshsizefilename );
} catch (netgen::NgException & ex) {
return error( COMPERR_BAD_PARMETERS, ex.What() );
}
}
TNodeToIDMap nodeToNetgenID;
nodeVec.clear();
ngMesh->AddFaceDescriptor( netgen::FaceDescriptor( faceID, solidID, solidID, 0 ));
// set local size according to size of existing segments
SMDS_ElemIteratorPtr iteratorElem = meshDS->elementsIterator(SMDSAbs_Edge);
while ( iteratorElem->more() ) // loop on elements on a geom face
{
const SMDS_MeshElement* elem = iteratorElem->next();
// Keeping only element that are in the element orientation file
bool isIn = elemOrientation.count(elem->GetID())==1;
if(!isIn)
continue;
bool isRev = elemOrientation[elem->GetID()];
std::cerr << isRev;
for ( int iN = 0; iN < 2; ++iN )
{
const SMDS_MeshNode* node = elem->GetNode( iN );
const int shapeID = node->getshapeId();
int& ngID = nodeToNetgenID.insert(TN2ID( node, invalid_ID )).first->second;
if ( ngID == invalid_ID )
{
ngID = ++Netgen_NbOfNodes;
Netgen_point [ 0 ] = node->X();
Netgen_point [ 1 ] = node->Y();
Netgen_point [ 2 ] = node->Z();
netgen::MeshPoint mp( netgen::Point<3> (node->X(), node->Y(), node->Z()) );
ngMesh->AddPoint ( mp, 1, netgen::EDGEPOINT );
}
Netgen_segment[ isRev ? 1-iN : iN ] = ngID;
}
// add segment
netgen::Segment seg;
seg[0] = Netgen_segment[0];
seg[1] = Netgen_segment[1];
seg.edgenr = ngMesh->GetNSeg() +1;
seg.si = faceID;
ngMesh->AddSegment(seg);
}
int nbNodes2 = ngMesh->GetNP();
int nseg = ngMesh->GetNSeg();
// insert old nodes into nodeVec
nodeVec.resize( nodeToNetgenID.size() + 1, 0 );
TNodeToIDMap::iterator n_id = nodeToNetgenID.begin();
for ( ; n_id != nodeToNetgenID.end(); ++n_id )
nodeVec[ n_id->second ] = n_id->first;
nodeToNetgenID.clear();
//if ( !isCommonLocalSize )
//limitSize( ngMesh, mparam.maxh * 0.8);
// -------------------------
// Generate surface mesh
// -------------------------
const int startWith = netgen::MESHCONST_MESHSURFACE;
const int endWith = toOptimize ? netgen::MESHCONST_OPTSURFACE : netgen::MESHCONST_MESHSURFACE;
SMESH_Comment str;
try {
OCC_CATCH_SIGNALS;
err = ngLib.GenerateMesh(occgeom, startWith, endWith, ngMesh);
if ( netgen::multithread.terminate )
return false;
if ( err )
str << "Error in netgen::OCCGenerateMesh() at " << netgen::multithread.task;
}
catch (Standard_Failure& ex)
{
err = 1;
str << "Exception in netgen::OCCGenerateMesh()"
<< " at " << netgen::multithread.task
<< ": " << ex.DynamicType()->Name();
if ( ex.GetMessageString() && strlen( ex.GetMessageString() ))
str << ": " << ex.GetMessageString();
}
catch (...) {
err = 1;
str << "Exception in netgen::OCCGenerateMesh()"
<< " at " << netgen::multithread.task;
}
if ( err )
{
if ( iLoop == LOC_SIZE )
{
std::cout << "Need second run" << std::endl;
/*netgen::mparam.minh = netgen::mparam.maxh;
netgen::mparam.maxh = 0;
for ( size_t iW = 0; iW < wires.size(); ++iW )
{
StdMeshers_FaceSidePtr wire = wires[ iW ];
const vector<UVPtStruct>& uvPtVec = wire->GetUVPtStruct();
for ( size_t iP = 1; iP < uvPtVec.size(); ++iP )
{
SMESH_TNodeXYZ p( uvPtVec[ iP ].node );
netgen::Point3d np( p.X(),p.Y(),p.Z());
double segLen = p.Distance( uvPtVec[ iP-1 ].node );
double size = ngMesh->GetH( np );
netgen::mparam.minh = Min( netgen::mparam.minh, size );
netgen::mparam.maxh = Max( netgen::mparam.maxh, segLen );
}
}
//cerr << "min " << mparam.minh << " max " << mparam.maxh << endl;
netgen::mparam.minh *= 0.9;
netgen::mparam.maxh *= 1.1;
*/
continue;
}
else
{
//faceErr.reset( new SMESH_ComputeError( COMPERR_ALGO_FAILED, str ));
}
}
// ----------------------------------------------------
// Fill the SMESHDS with the generated nodes and faces
// ----------------------------------------------------
if(output_mesh)
{
int nbNodes = ngMesh->GetNP();
int nbFaces = ngMesh->GetNSE();
std::cout << nbFaces << " " << nbNodes << std::endl;
int nbInputNodes = (int) nodeVec.size()-1;
nodeVec.resize( nbNodes+1, 0 );
// add nodes
for ( int ngID = nbInputNodes + 1; ngID <= nbNodes; ++ngID )
{
const netgen::MeshPoint& ngPoint = ngMesh->Point( ngID );
SMDS_MeshNode * node = meshDS->AddNode(ngPoint(0), ngPoint(1), ngPoint(2));
nodeVec[ ngID ] = node;
}
// create faces
int i,j;
for ( i = 1; i <= nbFaces ; ++i )
{
Ng_GetVolumeElement(ngLib.ngMesh(), i, Netgen_triangle);
helper.AddFace (nodeVec.at( Netgen_triangle[0] ),
nodeVec.at( Netgen_triangle[1] ),
nodeVec.at( Netgen_triangle[2] ));
}
} // output_mesh
break;
} // two attempts
//} // loop on FACEs
return true;
}

117
src/NETGENPlugin/NETGENPlugin_Runner.hxx
View File

@ -1,117 +0,0 @@
// Copyright (C) 2007-2021 CEA/DEN, EDF R&D, 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_Runner.hxx
// Author : Yoann AUDOUIN, EDF
// Module : NETGEN
//
#ifndef _NETGENPLUGIN_RUNNER_HXX_
#define _NETGENPLUGIN_RUNNER_HXX_
#include <string>
#include <iostream>
#include "NETGENPlugin_Defs.hxx"
#include "NETGENPlugin_Mesher.hxx"
#include "SMESH_Algo.hxx"
#include "Utils_SALOME_Exception.hxx"
class TopoDS_Shape;
class SMESH_Mesh;
class SMESH_Comment;
class netgen_params;
class NETGENPlugin_NetgenLibWrapper;
class SMDS_MeshNode;
// Netgen 2d functions
int netgen2dInternal(TopoDS_Shape &aShape,
SMESH_Mesh& aMesh,
netgen_params& aParams,
std::string new_element_file,
std::string element_orientation_file,
bool output_mesh);
int netgen2d(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);
// Netgen 3D functions
bool mycomputeFillNgMesh(
SMESH_Mesh& aMesh,
const TopoDS_Shape& aShape,
std::vector< const SMDS_MeshNode* > &nodeVec,
NETGENPlugin_NetgenLibWrapper &ngLib,
SMESH_MesherHelper &helper,
netgen_params &aParams,
std::string element_orientation_file,
int &Netgen_NbOfNodes);
bool mycomputePrepareParam(
SMESH_Mesh& aMesh,
NETGENPlugin_NetgenLibWrapper &ngLib,
netgen::OCCGeometry &occgeo,
SMESH_MesherHelper &helper,
netgen_params &aParams,
int &endWith);
bool mycomputeRunMesher(
netgen::OCCGeometry &occgeo,
std::vector< const SMDS_MeshNode* > &nodeVec,
netgen::Mesh* ngMesh,
NETGENPlugin_NetgenLibWrapper &ngLib,
int &startWith, int &endWith);
bool mycomputeFillNewElementFile(
std::vector< const SMDS_MeshNode* > &nodeVec,
NETGENPlugin_NetgenLibWrapper &ngLib,
std::string new_element_file,
int &Netgen_NbOfNodes);
bool mycomputeFillMesh(
std::vector< const SMDS_MeshNode* > &nodeVec,
NETGENPlugin_NetgenLibWrapper &ngLib,
SMESH_MesherHelper &helper,
int &Netgen_NbOfNodes);
int netgen3dInternal(TopoDS_Shape &aShape,
SMESH_Mesh& aMesh,
netgen_params& aParams,
std::string new_element_file,
std::string element_orientation_file,
bool output_mesh);
int netgen3d(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,
int nbThreads);
//TODO: Tmp function replace by real error handling
int error(int error_type, std::string msg);
int error(const SMESH_Comment& comment);
#endif

13
src/NETGENPlugin/NETGENPlugin_Runner_main.cxx
View File

@ -25,8 +25,7 @@
// Module : NETGEN
//
#include "NETGENPlugin_Runner.hxx"
#include "NETGENPlugin_NETGEN_3D_SA.hxx"
#include <stdio.h>
#include <string.h>
@ -84,7 +83,8 @@ int main(int argc, char *argv[]){
if (mesher=="NETGEN3D"){
auto begin = std::chrono::high_resolution_clock::now();
netgen3d(input_mesh_file,
NETGENPlugin_NETGEN_3D_SA myplugin;
myplugin.run(input_mesh_file,
shape_file,
hypo_file,
element_orientation_file,
@ -94,13 +94,6 @@ int main(int argc, char *argv[]){
auto end = std::chrono::high_resolution_clock::now();
auto elapsed = std::chrono::duration_cast<std::chrono::nanoseconds>(end - begin);
std::cout << "Time elapsed: " << elapsed.count()*1e-9 << std::endl;
} else if (mesher=="NETGEN2D"){
netgen2d(input_mesh_file,
shape_file,
hypo_file,
element_orientation_file,
new_element_file,
output_mesh_file);
} else {
std::cerr << "Unknown mesher:" << mesher << std::endl;
}

3
src/NETGENPlugin/NETGENPlugin_i.cxx
View File

@ -36,6 +36,7 @@
#include "NETGENPlugin_NETGEN_2D_ONLY_i.hxx"
#include "NETGENPlugin_NETGEN_2D_i.hxx"
#include "NETGENPlugin_NETGEN_3D_i.hxx"
#include "NETGENPlugin_NETGEN_3D_Remote_i.hxx"
#include "NETGENPlugin_SimpleHypothesis_2D_i.hxx"
#include "NETGENPlugin_SimpleHypothesis_3D_i.hxx"
@ -63,6 +64,8 @@ extern "C"
// Algorithms
if (strcmp(aHypName, "NETGEN_3D") == 0)
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") == 0)
aCreator = new NETGENPlugin_Creator_i<NETGENPlugin_NETGEN_2D_i>;
else if (strcmp(aHypName, "NETGEN_2D_ONLY") == 0)