NEMoSys  0.63.0
A modular, extensible resource with robust automated mesh generation, mesh quality analysis, adaptive mesh refinement, and data transfer between arbitrary meshes.
snappymeshGen Class Reference

snappymeshGen facilitates full-hexahedral/hex-dominent meshing of complex geometries with surface, patch, or region refinement capabilities in addition with adding layers to geometry. More...

Detailed Description

A background mesh is needed in order to perform meshing operation. Inputs can be given in form of surface files (.stl, .obj, etc.). Output mesh is written in OpenFOAM polyMesh format. snappymeshGen is based on snappyHexMesh utility of OpenFOAM with following features

  • Castallated Mesh
  • Snap
  • Add Layers

Definition at line 55 of file snappymeshGen.H.

Public Member Functions

 snappymeshGen ()
 snappymeshGen standard constructor More...
 
 snappymeshGen (snappymeshParams *params)
 snappymeshGen alternate constructor More...
 
 ~snappymeshGen ()
 snappymeshGen standard destructor More...
 
int createMeshFromSTL (const char *fname) override
 Creates mesh from input STL file. More...
 
vtkSmartPointer< vtkDataSet > getDataSet () const
 

Public Attributes

std::unique_ptr< NEM::MSH::geoMeshBasegmData
 

Protected Attributes

vtkSmartPointer< vtkDataSet > dataSet
 

Private Member Functions

void initialize ()
 Initializes OpenFOAM environment. More...
 
void createSnappyDict (const bool &write)
 Creates snappyHexMeshDict for meshing operation. More...
 

Private Attributes

bool defaults
 If enabled, generated mesh using default parameters. More...
 
snappymeshParamsparams_
 snappymeshParams object Parameters More...
 
std::unique_ptr< Foam::dictionary > controlDict_
 
std::unique_ptr< Foam::dictionary > fvSchemes_
 
std::unique_ptr< Foam::dictionary > fvSolution_
 
std::unique_ptr< Foam::dictionary > snappyMshDict_
 
std::unique_ptr< Foam::Time > runTime_
 
std::unique_ptr< Foam::fvMesh > fmesh_
 

Inherits meshGen.

Constructor & Destructor Documentation

◆ snappymeshGen() [1/2]

snappymeshGen::snappymeshGen ( )

Definition at line 49 of file snappymeshGen.C.

References defaults, initialize(), and params_.

50 {
51  // Default Meshing Parameters
52  params_ = new snappymeshParams();
53  defaults = true;
54 
55  // Initialization Tasks
56  initialize();
57 }
snappymeshParams contains all parameters essential for mesh generation using snappymeshGen class meth...
bool defaults
If enabled, generated mesh using default parameters.
Definition: snappymeshGen.H:89
snappymeshParams * params_
snappymeshParams object Parameters
Definition: snappymeshGen.H:93
void initialize()
Initializes OpenFOAM environment.
Definition: snappymeshGen.C:67

◆ snappymeshGen() [2/2]

snappymeshGen::snappymeshGen ( snappymeshParams params)
Parameters
paramssnappymeshParams object

Definition at line 60 of file snappymeshGen.C.

References initialize().

61  : defaults(false), params_(params) {
62  initialize(); // Foam Initialization
63 }
bool defaults
If enabled, generated mesh using default parameters.
Definition: snappymeshGen.H:89
snappymeshParams * params_
snappymeshParams object Parameters
Definition: snappymeshGen.H:93
void initialize()
Initializes OpenFOAM environment.
Definition: snappymeshGen.C:67

◆ ~snappymeshGen()

snappymeshGen::~snappymeshGen ( )

Definition at line 65 of file snappymeshGen.C.

65 {}

Member Function Documentation

◆ createMeshFromSTL()

int snappymeshGen::createMeshFromSTL ( const char *  fname)
overridevirtual

Implements meshGen.

Definition at line 91 of file snappymeshGen.C.

References fmesh_, meshGen::gmData, snappymeshParams::isPackMesh, params_, NEM::MSH::Read(), runTime_, and snappyMshDict_.

Referenced by NEM::DRV::SnappyMeshMeshGenDriver::execute().

91  {
92  auto objMsh = snappyMesh();
93 
94  bool writeMsh;
95  if (params_->isPackMesh)
96  writeMsh = true;
97  else
98  writeMsh = false;
99 
100  fmesh_ = std::unique_ptr<Foam::fvMesh>(
101  objMsh.generate(runTime_, snappyMshDict_, nullptr, writeMsh));
102 
103  // Create foamGeoMesh
104  if (writeMsh) {
105  gmData.reset();
106  gmData = std::unique_ptr<NEM::MSH::geoMeshBase>(NEM::MSH::Read(".foam"));
107  } else {
108  auto fgm_ = std::unique_ptr<NEM::MSH::foamGeoMesh>(
109  new NEM::MSH::foamGeoMesh(fmesh_.get(), ""));
110  gmData = std::unique_ptr<NEM::MSH::geoMeshBase>(
111  dynamic_cast<NEM::MSH::geoMeshBase *>(fgm_.get()));
112  }
113 
114  return 0;
115 }
std::unique_ptr< NEM::MSH::geoMeshBase > gmData
Definition: meshGen.H:54
geoMeshBase * Read(const std::string &fileName)
Read a mesh from file.
std::unique_ptr< Foam::Time > runTime_
bool isPackMesh
Boolean for packmesh.
std::unique_ptr< Foam::dictionary > snappyMshDict_
A concrete implementation of geoMeshBase representing a mesh in a fvMesh.
Definition: foamGeoMesh.H:52
std::unique_ptr< Foam::fvMesh > fmesh_
abstract class to specify geometry and mesh data
Definition: geoMeshBase.H:102
snappymeshParams * params_
snappymeshParams object Parameters
Definition: snappymeshGen.H:93

◆ createSnappyDict()

void snappymeshGen::createSnappyDict ( const bool &  write)
private

Definition at line 117 of file snappymeshGen.C.

References shmCastMeshControls::alwFreeZone, shmLayerControls::bufferCells, snappymeshParams::castMeshControls, shmCastMeshControls::castMeshGpLvl, shmCastMeshControls::cellsBetnLvls, shmMeshQualityControls::errReduction, shmSnapControls::explicitFeatureSnap, shmLayerControls::expRatio, shmCastMeshControls::featAngle, shmLayerControls::featAngle, shmLayerControls::finLThick, shmLayerControls::firstLyrThickness, shmCastMeshControls::ftrEdge, shmCastMeshControls::gapLvlInc, snappymeshParams::geoDef, snappymeshParams::geomFileName, shmCastMeshControls::geomRefs, shmSnapControls::implicitFeatureSnap, snappymeshParams::layerControls, shmLayerControls::layerVec, shmCastMeshControls::locMesh, shmMeshQualityControls::maxBndrySkew, shmMeshQualityControls::maxConc, shmLayerControls::maxFcTR, shmCastMeshControls::maxGCells, shmMeshQualityControls::maxIntSkew, shmCastMeshControls::maxLCells, shmMeshQualityControls::maxNonOrtho, shmLayerControls::maxThickTMR, snappymeshParams::mergeTol, shmMeshQualityControls::minArea, shmMeshQualityControls::minDet, shmMeshQualityControls::minFaceW, shmLayerControls::minMedAngl, shmCastMeshControls::minRefCells, shmMeshQualityControls::minTetQ, shmLayerControls::minThick, shmMeshQualityControls::minTriTwist, shmMeshQualityControls::minTwist, shmMeshQualityControls::minVol, shmMeshQualityControls::minVolRto, shmSnapControls::multiRegionFeatureSnap, shmSnapControls::nFeatureSnapIter, shmLayerControls::nGrow, shmLayerControls::nIter, shmLayerControls::nMedialAxisIter, shmLayerControls::nRelaxedIter, shmSnapControls::nRelaxIter, shmLayerControls::nSmoothDisplacement, shmSnapControls::nSmoothPatch, shmSnapControls::nSolveIter, params_, shmCastMeshControls::planarAngle, snappymeshParams::qualityControls, shmCastMeshControls::refSurfLvlMax, shmCastMeshControls::refSurfLvlMin, shmLayerControls::relaxIter, shmLayerControls::relSize, shmGeo::singleSolidPatch, shmLayerControls::slipFeatureAngle, shmMeshQualityControls::smoothScale, shmLayerControls::smthNorm, shmLayerControls::smthSurfNorm, shmLayerControls::smthThick, snappymeshParams::snapControls, snappyMshDict_, shmGeo::srchShape, shmGeo::stlPatchDefs, shmCastMeshControls::surfRefs, shmLayerControls::thickness, shmSnapControls::tolerance, snappymeshParams::withCastMesh, shmCastMeshControls::withCellZones, snappymeshParams::withLayers, shmGeo::withMultiPatches, and snappymeshParams::withSnap.

Referenced by initialize().

117  {
118  // header
119  std::string contText =
120  "\
121 /*--------------------------------*- C++ -*----------------------------------*\n\
122 | ========= | |\n\
123 | \\\\ / F ield | NEMoSys: snappyHexMesh interface |\n\
124 | \\\\ / O peration | |\n\
125 | \\\\ / A nd | |\n\
126 | \\\\/ M anipulation | |\n\
127 \\*---------------------------------------------------------------------------*/\n\
128 \n\
129 FoamFile\n\
130 {\n\
131  version 2.0;\n\
132  format ascii;\n\
133  class dictionary;\n\
134  location \"system\";\n\
135  object snappyHexMeshDict;\n\
136 }\n\n";
137 
138  contText =
139  contText +
140  "// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //\n\n";
141 
142  // initial booleans
143  if ((params_->withCastMesh)) {
144  contText += "\ncastellatedMesh true;\n";
145  } else {
146  contText += "\ncastellatedMesh false;\n";
147  }
148 
149  if ((params_->withSnap)) {
150  contText += "\nsnap true;\n";
151  } else {
152  contText += "\nsnap false;\n";
153  }
154 
155  if ((params_->withLayers)) {
156  contText += "\naddLayers true;\n";
157  } else {
158  contText += "\naddLayers false;\n";
159  }
160 
161  // Geometry declaration
162  contText += "\ngeometry\n{\n";
163  contText += "\t" + (params_->geomFileName);
164  contText += "\n\t{\n\t\ttype\ttriSurfaceMesh;\n";
165 
167  if (params_->geoDef.stlPatchDefs.empty()) {
168  std::cerr << "Error reading in custom patches from JSON input!"
169  << std::endl;
170  throw;
171  }
172  contText += "\t\tregions\n\t\t{\n";
173  for (auto pt = (params_->geoDef.stlPatchDefs).begin();
174  pt != (params_->geoDef.stlPatchDefs).end(); pt++) {
175  contText += "\t\t\t" + pt->STLPatchName + " \t{name " +
176  pt->snappyPatchName + ";}\n";
177  }
178  contText += "\t\t}\n";
179  } else {
180  contText =
181  contText + "\t\tname\t" + (params_->geoDef.singleSolidPatch) + ";\n";
182  }
183 
184  // Searchable shapes go here
185  if (!params_->geoDef.srchShape.empty()) {
186  for (const auto &shape : params_->geoDef.srchShape) {
187  contText += "\t\t" + (shape->patchName);
188  contText += "\n\t\t{";
189  contText += "\n\t\t\ttype " + (shape->getType()) + ";\n";
190 
191  if (auto box = std::dynamic_pointer_cast<shmSearchableBox>(shape)) {
192  contText += "\t\t\tmin (" + std::to_string(box->minBound.at(0)) + " " +
193  std::to_string(box->minBound.at(1)) + " " +
194  std::to_string(box->minBound.at(2)) + ");\n";
195  contText += "\t\t\tmax (" + std::to_string(box->maxBound.at(0)) + " " +
196  std::to_string(box->maxBound.at(1)) + " " +
197  std::to_string(box->maxBound.at(2)) + ");\n";
198  } else if (auto sphere =
199  std::dynamic_pointer_cast<shmSearchableSphere>(shape)) {
200  contText += "\t\t\tcentre (" + std::to_string(sphere->center.at(0)) +
201  " " + std::to_string(sphere->center.at(1)) + " " +
202  std::to_string(sphere->center.at(2)) + ");\n";
203  contText += "\t\t\tradius " + std::to_string(sphere->radius) + ";\n";
204  } else if (auto cyl =
205  std::dynamic_pointer_cast<shmSearchableCylinder>(shape)) {
206  contText += "\t\t\tpoint1 (" + std::to_string(cyl->axisPoint1.at(0)) +
207  " " + std::to_string(cyl->axisPoint1.at(1)) + " " +
208  std::to_string(cyl->axisPoint1.at(2)) + ");\n";
209  contText += "\t\t\tpoint2 (" + std::to_string(cyl->axisPoint2.at(0)) +
210  " " + std::to_string(cyl->axisPoint2.at(1)) + " " +
211  std::to_string(cyl->axisPoint2.at(2)) + ");\n";
212  contText += "\t\t\tradius " + std::to_string(cyl->radius) + ";\n";
213  }
214  contText += "\n\t\t}\n";
215  }
216  }
217  contText += "\t}\n";
218  contText += "};\n";
219  // End of geometry definition
220 
221  // Castellated Mesh Controls
222  contText += "\n\ncastellatedMeshControls\n{\n";
223  contText += "\tmaxLocalCells\t" +
224  std::to_string(params_->castMeshControls.maxLCells) + ";\n";
225  contText += "\tmaxGlobalCells\t" +
226  std::to_string(params_->castMeshControls.maxGCells) + ";\n";
227  contText += "\tnCellsBetweenLevels\t" +
228  std::to_string(params_->castMeshControls.cellsBetnLvls) + ";\n";
229  contText += "\tminRefinementCells\t" +
230  std::to_string(params_->castMeshControls.minRefCells) + ";\n";
231 
232  // Features using emesh file
233  contText += "\n\tfeatures\n\t(\n";
234  if (!params_->castMeshControls.ftrEdge.empty()) {
235  for (auto pt = (params_->castMeshControls.ftrEdge).begin();
236  pt != (params_->castMeshControls.ftrEdge).end(); pt++) {
237  contText += "\t\t{\n";
238  contText += "\t\t\tfile \"" + pt->fileName + "\";\n";
239  contText += "\t\t\tlevels ((" + std::to_string(pt->minLvl) + " " +
240  std::to_string(pt->maxLvl) + "));\n";
241  contText += "\t\t}\n";
242  }
243  }
244  contText += "\t);\n";
245 
246  contText += "\n\trefinementSurfaces\n\t{\n";
248  contText += "\t\t" + (params_->geomFileName) + "\n\t\t{";
249  else
250  contText =
251  contText + "\t\t" + (params_->geoDef.singleSolidPatch) + "\n\t\t{";
252  contText += "\n\t\t\tlevel (" +
253  std::to_string(params_->castMeshControls.refSurfLvlMin) + " " +
254  std::to_string(params_->castMeshControls.refSurfLvlMax) + ");\n";
255 
257  contText += "\t\t\tfaceZone\t" + (params_->geoDef.singleSolidPatch) + ";\n";
258  contText += "\t\t\tcellZone\t" + (params_->geoDef.singleSolidPatch) + ";\n";
259  contText += "\t\t\tcellZoneInside\tinside;\n";
260  }
261 
263  if ((!params_->castMeshControls.surfRefs.empty())) {
264  contText += "\n\t\t\tregions\n";
265  contText += "\t\t\t{";
266 
267  for (auto pt = (params_->castMeshControls.surfRefs).begin();
268  pt != (params_->castMeshControls.surfRefs).end(); pt++) {
269  if (pt->patchType == "NO") {
270  contText += "\n\t\t\t\t" + (pt->refPatchName) + "\t{level (" +
271  std::to_string((int)pt->minLvl) + " " +
272  std::to_string(pt->maxLvl) +
273  "); patchInfo { type patch; }}";
274  } else {
275  contText += "\n\t\t\t\t" + (pt->refPatchName) + "\t{level (" +
276  std::to_string((int)pt->minLvl) + " " +
277  std::to_string(pt->maxLvl) + "); patchInfo { type " +
278  pt->patchType + "; }}";
279  }
280  }
281 
282  contText += "\n\t\t\t}\n";
283  }
284  }
285 
286  contText += "\t\t}\n";
287 
288  contText += "\t\tgapLevelIncrement " +
289  std::to_string(params_->castMeshControls.castMeshGpLvl) + +";\n";
290 
291  contText += "\t}\n";
292 
293  contText += "\tresolveFeatureAngle\t" +
294  std::to_string(params_->castMeshControls.featAngle) + ";\n";
295 
296  contText += "\tgapLevelIncrement\t" +
297  std::to_string(params_->castMeshControls.gapLvlInc) + ";\n";
298 
299  contText += "\n\tplanarAngle " +
300  std::to_string(params_->castMeshControls.planarAngle) + ";\n";
301 
302  contText += "\n\trefinementRegions\n\t{\n";
303 
304  if ((!params_->castMeshControls.geomRefs.empty())) {
305  for (auto pt = (params_->castMeshControls.geomRefs).begin();
306  pt != (params_->castMeshControls.geomRefs).end(); pt++) {
307  contText += "\n\t\t" + (pt->patchName);
308  contText += "\n\t\t{";
309  contText += "\n\t\t\tmode " + (pt->mode) + ";\n";
310  contText += "\t\t\tlevels ((" + std::to_string(pt->minLvl) + " " +
311  std::to_string(pt->maxLvl) + "));\n";
312  contText += "\t\t}\n";
313  }
314  }
315 
316  contText += "\t}";
317 
318  contText += "\n\tlocationInMesh\t(" +
319  std::to_string(params_->castMeshControls.locMesh.at(0)) + " " +
320  std::to_string(params_->castMeshControls.locMesh.at(1)) + " " +
321  std::to_string(params_->castMeshControls.locMesh.at(2)) + ");\n";
323  contText += "\tallowFreeStandingZoneFaces\ttrue;\n";
325  contText += "\tallowFreeStandingZoneFaces\tfalse;\n";
326  contText += "}\n";
327 
328  // Snap Controls
329  contText += "\n\nsnapControls\n{\n";
330  contText += "\tnSmoothPatch\t" +
331  std::to_string(params_->snapControls.nSmoothPatch) + ";\n";
332  contText +=
333  "\ttolerance\t" + std::to_string(params_->snapControls.tolerance) + ";\n";
334  contText += "\tnSolveIter\t" +
335  std::to_string(params_->snapControls.nSolveIter) + ";\n";
336  contText += "\tnRelaxIter\t" +
337  std::to_string(params_->snapControls.nRelaxIter) + ";\n";
338 
339  // New Features
340  contText += "\tnFeatureSnapIter\t" +
341  std::to_string(params_->snapControls.nFeatureSnapIter) + ";\n";
342 
344  contText += "\timplicitFeatureSnap\t true;\n";
345  else
346  contText += "\timplicitFeatureSnap\t false;\n";
347 
349  contText += "\texplicitFeatureSnap\t true;\n";
350  else
351  contText += "\texplicitFeatureSnap\t false;\n";
352 
354  contText += "\tmultiRegionFeatureSnap\t true;\n";
355  else
356  contText += "\tmultiRegionFeatureSnap\t false;\n";
357 
358  contText += "}\n";
359 
360  // Layer Controls
361  contText += "\n\naddLayersControls\n{\n";
362  if ((params_->layerControls.relSize) == 1)
363  contText += "\trelativeSizes\ttrue;\n";
364  if ((params_->layerControls.relSize) == 0)
365  contText += "\trelativeSizes\tfalse;\n";
366 
367  // Layers
368  contText += "\tlayers\n\t{\n";
369  if (!params_->layerControls.layerVec.empty()) {
370  for (auto pt = (params_->layerControls.layerVec).begin();
371  pt != (params_->layerControls.layerVec).end(); pt++) {
372  contText += "\t\t" + pt->patchName + "\n\t\t{\n";
373  contText +=
374  "\t\t\tnSurfaceLayers " + std::to_string(pt->nSurfaceLayers) + ";\n";
375  contText +=
376  "\t\t\texpansionRatio " + std::to_string(pt->expansionRatio) + ";\n";
377  contText += "\t\t\tfinalLayerThickness " +
378  std::to_string(pt->finalLayerThickness) + ";\n";
379  if (pt->firstLyrThickness > 0)
380  contText += "\t\t\tfirstLayerThickness " +
381  std::to_string(pt->firstLyrThickness) + ";\n";
382  if (pt->thickness > 0)
383  contText += "\t\t\tthickness " + std::to_string(pt->thickness) + ";\n";
384  contText +=
385  "\t\t\tminThickness " + std::to_string(pt->minThickness) + ";\n";
386  contText += "\t\t}\n";
387  }
388  }
389  contText += "\t}\n";
390 
391  contText += "\texpansionRatio\t" +
392  std::to_string(params_->layerControls.expRatio) + ";\n";
393  contText += "\tfinalLayerThickness\t" +
394  std::to_string(params_->layerControls.finLThick) + ";\n";
395  contText += "\tminThickness\t" +
396  std::to_string(params_->layerControls.minThick) + ";\n";
398  contText += "\tfirstLayerThickness\t" +
399  std::to_string(params_->layerControls.firstLyrThickness) +
400  ";\n";
402  contText += "\tthickness\t" +
403  std::to_string(params_->layerControls.thickness) + ";\n";
404  contText +=
405  "\tnGrow\t" + std::to_string(params_->layerControls.nGrow) + ";\n";
406  contText += "\tfeatureAngle\t" +
407  std::to_string(params_->layerControls.featAngle) + ";\n";
408  contText += "\tslipFeatureAngle\t" +
409  std::to_string(params_->layerControls.slipFeatureAngle) + ";\n";
410  contText += "\tnRelaxIter\t" +
411  std::to_string(params_->layerControls.relaxIter) + ";\n";
412  contText += "\tnSmoothSurfaceNormals\t" +
413  std::to_string(params_->layerControls.smthSurfNorm) + ";\n";
414  contText += "\tnSmoothNormals\t" +
415  std::to_string(params_->layerControls.smthNorm) + ";\n";
416  contText += "\tnSmoothThickness\t" +
417  std::to_string(params_->layerControls.smthThick) + ";\n";
418  contText += "\tmaxFaceThicknessRatio\t" +
419  std::to_string(params_->layerControls.maxFcTR) + ";\n";
420  contText += "\tmaxThicknessToMedialRatio\t" +
421  std::to_string(params_->layerControls.maxThickTMR) + ";\n";
422  contText += "\tminMedialAxisAngle\t" +
423  std::to_string(params_->layerControls.minMedAngl) + ";\n";
424  contText += "\tnBufferCellsNoExtrude\t" +
425  std::to_string(params_->layerControls.bufferCells) + ";\n";
426  contText +=
427  "\tnLayerIter\t" + std::to_string(params_->layerControls.nIter) + ";\n";
428  contText += "\tnRelaxedIter\t" +
429  std::to_string(params_->layerControls.nRelaxedIter) + ";\n";
431  contText += "\tnMedialAxisIter\t" +
432  std::to_string(params_->layerControls.nMedialAxisIter) + ";\n";
434  contText += "\tnSmoothDisplacement\t" +
435  std::to_string(params_->layerControls.nSmoothDisplacement) +
436  ";\n";
437  contText += "}\n";
438 
439  // Mesh Quality Controls
440  std::ostringstream minimumVol;
441  minimumVol << params_->qualityControls.minVol;
442  std::ostringstream minimumTetQuality;
443  minimumTetQuality << params_->qualityControls.minTetQ;
444 
445  contText += "\n\nmeshQualityControls\n{\n";
446  contText += "\tmaxNonOrtho\t" +
447  std::to_string(params_->qualityControls.maxNonOrtho) + ";\n";
448  contText += "\tmaxBoundarySkewness\t" +
449  std::to_string(params_->qualityControls.maxBndrySkew) + ";\n";
450  contText += "\tmaxInternalSkewness\t" +
451  std::to_string(params_->qualityControls.maxIntSkew) + ";\n";
452  contText += "\tmaxConcave\t" +
453  std::to_string(params_->qualityControls.maxConc) + ";\n";
454  contText += "\tminVol\t" + minimumVol.str() + ";\n";
455  contText += "\tminTetQuality\t" + minimumTetQuality.str() + ";\n";
456  contText +=
457  "\tminArea\t" + std::to_string(params_->qualityControls.minArea) + ";\n";
458  contText += "\tminTwist\t" +
459  std::to_string(params_->qualityControls.minTwist) + ";\n";
460  contText += "\tminFaceWeight\t" +
461  std::to_string(params_->qualityControls.minFaceW) + ";\n";
462  contText += "\tminVolRatio\t" +
463  std::to_string(params_->qualityControls.minVolRto) + ";\n";
464  contText += "\tminDeterminant\t" +
465  std::to_string(params_->qualityControls.minDet) + ";\n";
466  contText += "\tminTriangleTwist\t" +
467  std::to_string(params_->qualityControls.minTriTwist) + ";\n";
468  contText += "\tnSmoothScale\t" +
469  std::to_string(params_->qualityControls.smoothScale) + ";\n";
470  contText += "\terrorReduction\t" +
471  std::to_string(params_->qualityControls.errReduction) + ";\n";
472  contText += "\trelaxed\n\t{\n";
473  contText += "\t\tmaxNonOrtho \t " +
474  std::to_string(params_->qualityControls.maxNonOrtho) + ";\n\t}\n";
475  contText += "}\n";
476 
477  contText +=
478  "\nmergeTolerance\t" + std::to_string(params_->mergeTol) + ";\n\n";
479 
480  contText += "writeFlags\n(\n";
481  contText += "\tscalarLevels\n";
482  contText += "\tlayerSets\n";
483  contText += "\tlayerFields\n);\n";
484 
485  contText =
486  contText +
487  "// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //";
488 
489  // Keep this dictionary in memory
490  Foam::dictionary tmptmpDuc =
491  new Foam::dictionary(Foam::IStringStream(contText)(), true);
492  snappyMshDict_ = std::unique_ptr<Foam::dictionary>(
493  new Foam::dictionary("snappyHexMeshDict"));
494  snappyMshDict_->merge(tmptmpDuc);
495 
496  // Write in traditional way due to formatting issues
497  if (write) {
498  // creating a base system directory
499  const char dir_path[] = "./system";
500  boost::filesystem::path dir(dir_path);
501  try {
502  boost::filesystem::create_directory(dir);
503  } catch (boost::filesystem::filesystem_error &e) {
504  std::cerr << "Problem in creating system directory for the snappyHexMesh"
505  << "\n";
506  std::cerr << e.what() << std::endl;
507  throw;
508  }
509 
510  // creating mesh dictionary file
511  std::ofstream contDict;
512  contDict.open(std::string(dir_path) + "/snappyHexMeshDict");
513  contDict << contText;
514  contDict.close();
515  }
516 }
int maxGCells
max local cells allowed (on 1 processor)
std::string singleSolidPatch
patch name for input surface file for when not _withMultiPatches
std::vector< shmSurfRefine > surfRefs
Vector for shmSurfRefine struct.
int nMedialAxisIter
Limits # of steps walking away from the surface.
int relaxIter
number of relaxation interations for layer addition
int cellsBetnLvls
number of cells between levels
double maxConc
Max Concativity.
shmMeshQualityControls qualityControls
Mesh Quality Controls.
double minVolRto
Minimum Volume Ratio.
double featAngle
resolve feature angle
double expRatio
Expansion Ratio for layer addition.
int nIter
Number of layer interations for layer addition.
int nFeatureSnapIter
Feature snapping.
shmGeo geoDef
Geometry definition.
int planarAngle
Planar angle castallated mesh.
bool withSnap
Enables surface snapping option for snappymeshGen.
int maxNonOrtho
Maximum non-orthogonality.
double minThick
Minimum Thickness for layer addition.
bool multiRegionFeatureSnap
Multi-region feature snapping boolean.
int gapLvlInc
Gap level increment.
std::vector< shmSTLDefinition > stlPatchDefs
Vector for STL patches, if _withMultiPatches.
double firstLyrThickness
First Layer Thickness.
shmSnapControls snapControls
Snapping controls.
int nSmoothPatch
Number of smoothing patches during snapping procedure.
double maxBndrySkew
Max Boundary Skewness.
double tolerance
Snapping tolerance at shared interface (affects conformality of mesh)
double maxThickTMR
Maximum thickness to medial ratio for layer addition.
std::vector< shmFeatureEdgeRef > ftrEdge
Vector for feature edge refinement.
std::string geomFileName
Input geometry STL name.
std::vector< shmLayers > layerVec
Vector for Layers.
bool implicitFeatureSnap
Implicit feature snaping boolean.
double minFaceW
Minimum Face Weight.
double minMedAngl
Minimum medial axis angle for layer addition.
bool relSize
Enables relative sizes option during layer addition.
std::vector< std::shared_ptr< shmSearchableShape > > srchShape
Vector for searchable shapes.
double minTwist
Minimum Twist.
int smthNorm
Number of smooth normals for layer addition.
double thickness
Thickness.
bool withMultiPatches
Boolean for if user wants whole STL as one solid of has different patches under STL.
double maxIntSkew
Max Internal Skewness.
double minVol
Minimum Cell Volume.
int nRelaxedIter
Relaxed iteration.
int smthThick
Number of smooth thickness for layer addition.
int bufferCells
Number of buffer cells no extrude for layer addition.
bool alwFreeZone
allows free standing zones if enabled
double minArea
Minimum Area.
int smoothScale
nSmoothScale
bool withLayers
Enables add layers option for snappymeshGen.
bool withCellZones
Enables defining cellzones and facezones.
std::unique_ptr< Foam::dictionary > snappyMshDict_
int minRefCells
minimum refinement cells
shmCastMeshControls castMeshControls
Castellated mesh controls.
bool explicitFeatureSnap
Explicit feature snapping boolean.
int maxLCells
max global cells allowed in mesh
double errReduction
Error Reduction.
double minTriTwist
Minimum Triangle Twist.
std::vector< shmRegionRefine > geomRefs
Vector for shmRegionRefine struct.
double mergeTol
merge tolerance for mesh
int nSolveIter
Maximum iterations during snapping procedure before ending the process.
shmLayerControls layerControls
Layer controls.
double maxFcTR
Maximum face thickness ratio for layer addition.
std::array< double, 3 > locMesh
location in Mesh (cells are kept in regions accessinble from this location)
double minDet
Minimum Determinant.
int nRelaxIter
Maximum number of relaxation iterations for snapping procedure.
double featAngle
Feature Angle for layer addition.
double minTetQ
Minimum Tet Quality.
int castMeshGpLvl
Gap level increment.
int smthSurfNorm
Number of smooth surface normals for layer addition.
int nGrow
Growth rate of successive layers.
int refSurfLvlMin
minimum surface refinement
double finLThick
Final Layer Thickness for layer addition.
int slipFeatureAngle
Slip feature angles.
bool withCastMesh
Enables castellated mesh option for snappymeshGen.
snappymeshParams * params_
snappymeshParams object Parameters
Definition: snappymeshGen.H:93
int nSmoothDisplacement
Smooth displacement after medial axis determination.
int refSurfLvlMax
maximum surface refinement

◆ getDataSet()

vtkSmartPointer<vtkDataSet> meshGen::getDataSet ( ) const
inlineinherited

Definition at line 53 of file meshGen.H.

Referenced by NEM::DRV::MeshGenDriver::MeshGenDriver().

53 { return dataSet; }
vtkSmartPointer< vtkDataSet > dataSet
Definition: meshGen.H:57

◆ initialize()

void snappymeshGen::initialize ( )
private

Definition at line 67 of file snappymeshGen.C.

References controlDict_, createSnappyDict(), fvSchemes_, fvSolution_, snappymeshParams::isPackMesh, params_, and runTime_.

Referenced by snappymeshGen().

67  {
68  // create dictionaries needed
69  bool writeDicts;
70  if (params_->isPackMesh)
71  writeDicts = true;
72  else
73  writeDicts = false;
74 
75  std::unique_ptr<getDicts> initFoam;
76  initFoam = std::unique_ptr<getDicts>(new getDicts());
77  controlDict_ = initFoam->createControlDict(writeDicts);
78  fvSchemes_ = initFoam->createFvSchemes(writeDicts);
79  fvSolution_ = initFoam->createFvSolution(writeDicts);
80 
81  // create dictionaries needed in memory
82  createSnappyDict(writeDicts);
83 
84  runTime_ =
85  std::unique_ptr<Foam::Time>(new Foam::Time(controlDict_.get(), ".", "."));
86 
87  //- 2d cartesian mesher cannot be run in parallel
88  Foam::argList::noParallel();
89 }
std::unique_ptr< Foam::Time > runTime_
std::unique_ptr< Foam::dictionary > fvSchemes_
Definition: snappymeshGen.H:98
bool isPackMesh
Boolean for packmesh.
std::unique_ptr< Foam::dictionary > fvSolution_
Definition: snappymeshGen.H:99
void createSnappyDict(const bool &write)
Creates snappyHexMeshDict for meshing operation.
snappymeshParams * params_
snappymeshParams object Parameters
Definition: snappymeshGen.H:93
std::unique_ptr< Foam::dictionary > controlDict_
Definition: snappymeshGen.H:97

Member Data Documentation

◆ controlDict_

std::unique_ptr<Foam::dictionary> snappymeshGen::controlDict_
private

Definition at line 97 of file snappymeshGen.H.

Referenced by initialize().

◆ dataSet

vtkSmartPointer<vtkDataSet> meshGen::dataSet
protectedinherited

Definition at line 57 of file meshGen.H.

◆ defaults

bool snappymeshGen::defaults
private

Definition at line 89 of file snappymeshGen.H.

Referenced by snappymeshGen().

◆ fmesh_

std::unique_ptr<Foam::fvMesh> snappymeshGen::fmesh_
private

Definition at line 102 of file snappymeshGen.H.

Referenced by createMeshFromSTL().

◆ fvSchemes_

std::unique_ptr<Foam::dictionary> snappymeshGen::fvSchemes_
private

Definition at line 98 of file snappymeshGen.H.

Referenced by initialize().

◆ fvSolution_

std::unique_ptr<Foam::dictionary> snappymeshGen::fvSolution_
private

Definition at line 99 of file snappymeshGen.H.

Referenced by initialize().

◆ gmData

◆ params_

snappymeshParams* snappymeshGen::params_
private

Definition at line 93 of file snappymeshGen.H.

Referenced by createMeshFromSTL(), createSnappyDict(), initialize(), and snappymeshGen().

◆ runTime_

std::unique_ptr<Foam::Time> snappymeshGen::runTime_
private

Definition at line 101 of file snappymeshGen.H.

Referenced by createMeshFromSTL(), and initialize().

◆ snappyMshDict_

std::unique_ptr<Foam::dictionary> snappymeshGen::snappyMshDict_
private

Definition at line 100 of file snappymeshGen.H.

Referenced by createMeshFromSTL(), and createSnappyDict().


The documentation for this class was generated from the following files: