snappyHexMesh.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2019 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 Application
25  snappyHexMesh
26 
27 Description
28  Automatic split hex mesher. Refines and snaps to surface.
29 
30 \*---------------------------------------------------------------------------*/
31 
32 #include "argList.H"
33 #include "Time.H"
34 #include "fvMesh.H"
35 #include "snappyRefineDriver.H"
36 #include "snappySnapDriver.H"
37 #include "snappyLayerDriver.H"
38 #include "searchableSurfaces.H"
39 #include "refinementSurfaces.H"
40 #include "refinementFeatures.H"
41 #include "shellSurfaces.H"
42 #include "decompositionMethod.H"
43 #include "noDecomp.H"
44 #include "fvMeshDistribute.H"
45 #include "wallPolyPatch.H"
46 #include "refinementParameters.H"
47 #include "snapParameters.H"
48 #include "layerParameters.H"
49 #include "vtkSetWriter.H"
50 #include "faceSet.H"
51 #include "motionSmoother.H"
52 #include "polyTopoChange.H"
53 #include "cellModeller.H"
55 #include "surfZoneIdentifierList.H"
56 #include "UnsortedMeshedSurface.H"
57 #include "MeshedSurface.H"
58 #include "globalIndex.H"
59 #include "IOmanip.H"
60 #include "fvMeshTools.H"
61 
62 using namespace Foam;
63 
64 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
65 
66 // Convert size (as fraction of defaultCellSize) to refinement level
67 label sizeCoeffToRefinement
68 (
69  const scalar level0Coeff, // ratio of hex cell size v.s. defaultCellSize
70  const scalar sizeCoeff
71 )
72 {
73  return round(::log(level0Coeff/sizeCoeff)/::log(2));
74 }
75 
76 
77 autoPtr<refinementSurfaces> createRefinementSurfaces
78 (
79  const searchableSurfaces& allGeometry,
80  const dictionary& surfacesDict,
81  const dictionary& shapeControlDict,
82  const label gapLevelIncrement,
83  const scalar level0Coeff
84 )
85 {
86  autoPtr<refinementSurfaces> surfacePtr;
87 
88  // Count number of surfaces.
89  label surfI = 0;
90  forAll(allGeometry.names(), geomI)
91  {
92  const word& geomName = allGeometry.names()[geomI];
93 
94  if (surfacesDict.found(geomName))
95  {
96  surfI++;
97  }
98  }
99 
100  labelList surfaces(surfI);
101  wordList names(surfI);
102  PtrList<surfaceZonesInfo> surfZones(surfI);
103 
104  labelList regionOffset(surfI);
105 
106  labelList globalMinLevel(surfI, 0);
107  labelList globalMaxLevel(surfI, 0);
108  labelList globalLevelIncr(surfI, 0);
109  PtrList<dictionary> globalPatchInfo(surfI);
110  List<Map<label>> regionMinLevel(surfI);
111  List<Map<label>> regionMaxLevel(surfI);
112  List<Map<label>> regionLevelIncr(surfI);
113  List<Map<scalar>> regionAngle(surfI);
114  List<Map<autoPtr<dictionary>>> regionPatchInfo(surfI);
115 
116  HashSet<word> unmatchedKeys(surfacesDict.toc());
117 
118  surfI = 0;
119  forAll(allGeometry.names(), geomI)
120  {
121  const word& geomName = allGeometry.names()[geomI];
122 
123  const entry* ePtr = surfacesDict.lookupEntryPtr(geomName, false, true);
124 
125  if (ePtr)
126  {
127  const dictionary& shapeDict = ePtr->dict();
128  unmatchedKeys.erase(ePtr->keyword());
129 
130  names[surfI] = geomName;
131  surfaces[surfI] = geomI;
132 
133  const searchableSurface& surface = allGeometry[geomI];
134 
135  // Find the index in shapeControlDict
136  // Invert surfaceCellSize to get the refinementLevel
137 
138  const word scsFuncName =
139  shapeDict.lookup("surfaceCellSizeFunction");
140  const dictionary& scsDict =
141  shapeDict.optionalSubDict(scsFuncName + "Coeffs");
142 
143  const scalar surfaceCellSize =
144  scsDict.lookup<scalar>("surfaceCellSizeCoeff");
145 
146  const label refLevel = sizeCoeffToRefinement
147  (
148  level0Coeff,
149  surfaceCellSize
150  );
151 
152  globalMinLevel[surfI] = refLevel;
153  globalMaxLevel[surfI] = refLevel;
154  globalLevelIncr[surfI] = gapLevelIncrement;
155 
156  // Surface zones
157  surfZones.set(surfI, new surfaceZonesInfo(surface, shapeDict));
158 
159 
160  // Global perpendicular angle
161  if (shapeDict.found("patchInfo"))
162  {
163  globalPatchInfo.set
164  (
165  surfI,
166  shapeDict.subDict("patchInfo").clone()
167  );
168  }
169 
170 
171  // Per region override of patchInfo
172 
173  if (shapeDict.found("regions"))
174  {
175  const dictionary& regionsDict = shapeDict.subDict("regions");
176  const wordList& regionNames =
177  allGeometry[surfaces[surfI]].regions();
178 
179  forAll(regionNames, regionI)
180  {
181  if (regionsDict.found(regionNames[regionI]))
182  {
183  // Get the dictionary for region
184  const dictionary& regionDict = regionsDict.subDict
185  (
186  regionNames[regionI]
187  );
188 
189  if (regionDict.found("patchInfo"))
190  {
191  regionPatchInfo[surfI].insert
192  (
193  regionI,
194  regionDict.subDict("patchInfo").clone()
195  );
196  }
197  }
198  }
199  }
200 
201  // Per region override of cellSize
202  if (shapeDict.found("regions"))
203  {
204  const dictionary& shapeControlRegionsDict =
205  shapeDict.subDict("regions");
206  const wordList& regionNames =
207  allGeometry[surfaces[surfI]].regions();
208 
209  forAll(regionNames, regionI)
210  {
211  if (shapeControlRegionsDict.found(regionNames[regionI]))
212  {
213  const dictionary& shapeControlRegionDict =
214  shapeControlRegionsDict.subDict
215  (
216  regionNames[regionI]
217  );
218 
219  const word scsFuncName =
220  shapeControlRegionDict.lookup
221  (
222  "surfaceCellSizeFunction"
223  );
224  const dictionary& scsDict =
225  shapeControlRegionDict.subDict
226  (
227  scsFuncName + "Coeffs"
228  );
229 
230  const scalar surfaceCellSize =
231  scsDict.lookup<scalar>("surfaceCellSizeCoeff");
232 
233  const label refLevel = sizeCoeffToRefinement
234  (
235  level0Coeff,
236  surfaceCellSize
237  );
238 
239  regionMinLevel[surfI].insert(regionI, refLevel);
240  regionMaxLevel[surfI].insert(regionI, refLevel);
241  regionLevelIncr[surfI].insert(regionI, 0);
242  }
243  }
244  }
245 
246  surfI++;
247  }
248  }
249 
250  // Calculate local to global region offset
251  label nRegions = 0;
252 
253  forAll(surfaces, surfI)
254  {
255  regionOffset[surfI] = nRegions;
256  nRegions += allGeometry[surfaces[surfI]].regions().size();
257  }
258 
259  // Rework surface specific information into information per global region
260  labelList minLevel(nRegions, 0);
261  labelList maxLevel(nRegions, 0);
262  labelList gapLevel(nRegions, -1);
263  PtrList<dictionary> patchInfo(nRegions);
264 
265  forAll(globalMinLevel, surfI)
266  {
267  label nRegions = allGeometry[surfaces[surfI]].regions().size();
268 
269  // Initialise to global (i.e. per surface)
270  for (label i = 0; i < nRegions; i++)
271  {
272  label globalRegionI = regionOffset[surfI] + i;
273  minLevel[globalRegionI] = globalMinLevel[surfI];
274  maxLevel[globalRegionI] = globalMaxLevel[surfI];
275  gapLevel[globalRegionI] =
276  maxLevel[globalRegionI]
277  + globalLevelIncr[surfI];
278 
279  if (globalPatchInfo.set(surfI))
280  {
281  patchInfo.set
282  (
283  globalRegionI,
284  globalPatchInfo[surfI].clone()
285  );
286  }
287  }
288 
289  // Overwrite with region specific information
290  forAllConstIter(Map<label>, regionMinLevel[surfI], iter)
291  {
292  label globalRegionI = regionOffset[surfI] + iter.key();
293 
294  minLevel[globalRegionI] = iter();
295  maxLevel[globalRegionI] = regionMaxLevel[surfI][iter.key()];
296  gapLevel[globalRegionI] =
297  maxLevel[globalRegionI]
298  + regionLevelIncr[surfI][iter.key()];
299  }
300 
301  const Map<autoPtr<dictionary>>& localInfo = regionPatchInfo[surfI];
302  forAllConstIter(Map<autoPtr<dictionary>>, localInfo, iter)
303  {
304  label globalRegionI = regionOffset[surfI] + iter.key();
305  patchInfo.set(globalRegionI, iter()().clone());
306  }
307  }
308 
309  surfacePtr.set
310  (
312  (
313  allGeometry,
314  surfaces,
315  names,
316  surfZones,
317  regionOffset,
318  minLevel,
319  maxLevel,
320  gapLevel,
321  scalarField(nRegions, -great), // perpendicularAngle,
322  patchInfo
323  )
324  );
325 
326 
327  const refinementSurfaces& rf = surfacePtr();
328 
329  // Determine maximum region name length
330  label maxLen = 0;
331  forAll(rf.surfaces(), surfI)
332  {
333  label geomI = rf.surfaces()[surfI];
334  const wordList& regionNames = allGeometry.regionNames()[geomI];
335  forAll(regionNames, regionI)
336  {
337  maxLen = Foam::max(maxLen, label(regionNames[regionI].size()));
338  }
339  }
340 
341 
342  Info<< setw(maxLen) << "Region"
343  << setw(10) << "Min Level"
344  << setw(10) << "Max Level"
345  << setw(10) << "Gap Level" << nl
346  << setw(maxLen) << "------"
347  << setw(10) << "---------"
348  << setw(10) << "---------"
349  << setw(10) << "---------" << endl;
350 
351  forAll(rf.surfaces(), surfI)
352  {
353  label geomI = rf.surfaces()[surfI];
354 
355  Info<< rf.names()[surfI] << ':' << nl;
356 
357  const wordList& regionNames = allGeometry.regionNames()[geomI];
358 
359  forAll(regionNames, regionI)
360  {
361  label globalI = rf.globalRegion(surfI, regionI);
362 
363  Info<< setw(maxLen) << regionNames[regionI]
364  << setw(10) << rf.minLevel()[globalI]
365  << setw(10) << rf.maxLevel()[globalI]
366  << setw(10) << rf.gapLevel()[globalI] << endl;
367  }
368  }
369 
370 
371  return surfacePtr;
372 }
373 
374 
375 void extractSurface
376 (
377  const polyMesh& mesh,
378  const Time& runTime,
379  const labelHashSet& includePatches,
380  const fileName& outFileName
381 )
382 {
383  const polyBoundaryMesh& bMesh = mesh.boundaryMesh();
384 
385  // Collect sizes. Hash on names to handle local-only patches (e.g.
386  // processor patches)
387  HashTable<label> patchSize(1000);
388  label nFaces = 0;
389  forAllConstIter(labelHashSet, includePatches, iter)
390  {
391  const polyPatch& pp = bMesh[iter.key()];
392  patchSize.insert(pp.name(), pp.size());
393  nFaces += pp.size();
394  }
396 
397 
398  // Allocate zone/patch for all patches
399  HashTable<label> compactZoneID(1000);
400  forAllConstIter(HashTable<label>, patchSize, iter)
401  {
402  label sz = compactZoneID.size();
403  compactZoneID.insert(iter.key(), sz);
404  }
405  Pstream::mapCombineScatter(compactZoneID);
406 
407 
408  // Rework HashTable into labelList just for speed of conversion
409  labelList patchToCompactZone(bMesh.size(), -1);
410  forAllConstIter(HashTable<label>, compactZoneID, iter)
411  {
412  label patchi = bMesh.findPatchID(iter.key());
413  if (patchi != -1)
414  {
415  patchToCompactZone[patchi] = iter();
416  }
417  }
418 
419  // Collect faces on zones
420  DynamicList<label> faceLabels(nFaces);
421  DynamicList<label> compactZones(nFaces);
422  forAllConstIter(labelHashSet, includePatches, iter)
423  {
424  const polyPatch& pp = bMesh[iter.key()];
425  forAll(pp, i)
426  {
427  faceLabels.append(pp.start()+i);
428  compactZones.append(patchToCompactZone[pp.index()]);
429  }
430  }
431 
432  // Addressing engine for all faces
433  uindirectPrimitivePatch allBoundary
434  (
435  UIndirectList<face>(mesh.faces(), faceLabels),
436  mesh.points()
437  );
438 
439 
440  // Find correspondence to master points
441  labelList pointToGlobal;
442  labelList uniqueMeshPoints;
443  autoPtr<globalIndex> globalNumbers = mesh.globalData().mergePoints
444  (
445  allBoundary.meshPoints(),
446  allBoundary.meshPointMap(),
447  pointToGlobal,
448  uniqueMeshPoints
449  );
450 
451  // Gather all unique points on master
452  List<pointField> gatheredPoints(Pstream::nProcs());
453  gatheredPoints[Pstream::myProcNo()] = pointField
454  (
455  mesh.points(),
456  uniqueMeshPoints
457  );
458  Pstream::gatherList(gatheredPoints);
459 
460  // Gather all faces
461  List<faceList> gatheredFaces(Pstream::nProcs());
462  gatheredFaces[Pstream::myProcNo()] = allBoundary.localFaces();
463  forAll(gatheredFaces[Pstream::myProcNo()], i)
464  {
465  inplaceRenumber(pointToGlobal, gatheredFaces[Pstream::myProcNo()][i]);
466  }
467  Pstream::gatherList(gatheredFaces);
468 
469  // Gather all ZoneIDs
470  List<labelList> gatheredZones(Pstream::nProcs());
471  gatheredZones[Pstream::myProcNo()] = move(compactZones);
472  Pstream::gatherList(gatheredZones);
473 
474  // On master combine all points, faces, zones
475  if (Pstream::master())
476  {
477  pointField allPoints = ListListOps::combine<pointField>
478  (
479  gatheredPoints,
481  );
482  gatheredPoints.clear();
483 
484  faceList allFaces = ListListOps::combine<faceList>
485  (
486  gatheredFaces,
488  );
489  gatheredFaces.clear();
490 
491  labelList allZones = ListListOps::combine<labelList>
492  (
493  gatheredZones,
495  );
496  gatheredZones.clear();
497 
498 
499  // Zones
500  surfZoneIdentifierList surfZones(compactZoneID.size());
501  forAllConstIter(HashTable<label>, compactZoneID, iter)
502  {
503  surfZones[iter()] = surfZoneIdentifier(iter.key(), iter());
504  Info<< "surfZone " << iter() << " : " << surfZones[iter()].name()
505  << endl;
506  }
507 
508  UnsortedMeshedSurface<face> unsortedFace
509  (
510  move(allPoints),
511  move(allFaces),
512  move(allZones),
513  move(surfZones)
514  );
515 
516 
517  MeshedSurface<face> sortedFace(unsortedFace);
518 
519  fileName globalCasePath
520  (
521  runTime.processorCase()
522  ? runTime.path()/".."/outFileName
523  : runTime.path()/outFileName
524  );
525  globalCasePath.clean();
526 
527  Info<< "Writing merged surface to " << globalCasePath << endl;
528 
529  sortedFace.write(globalCasePath);
530  }
531 }
532 
533 
534 // Check writing tolerance before doing any serious work
535 scalar getMergeDistance(const polyMesh& mesh, const scalar mergeTol)
536 {
537  const boundBox& meshBb = mesh.bounds();
538  scalar mergeDist = mergeTol * meshBb.mag();
539 
540  Info<< nl
541  << "Overall mesh bounding box : " << meshBb << nl
542  << "Relative tolerance : " << mergeTol << nl
543  << "Absolute matching distance : " << mergeDist << nl
544  << endl;
545 
546  // check writing tolerance
547  if (mesh.time().writeFormat() == IOstream::ASCII)
548  {
549  const scalar writeTol = std::pow
550  (
551  scalar(10.0),
552  -scalar(IOstream::defaultPrecision())
553  );
554 
555  if (mergeTol < writeTol)
556  {
558  << "Your current settings specify ASCII writing with "
559  << IOstream::defaultPrecision() << " digits precision." << nl
560  << "Your merging tolerance (" << mergeTol
561  << ") is finer than this." << nl
562  << "Change to binary writeFormat, "
563  << "or increase the writePrecision" << endl
564  << "or adjust the merge tolerance (mergeTol)."
565  << exit(FatalError);
566  }
567  }
568 
569  return mergeDist;
570 }
571 
572 
573 void removeZeroSizedPatches(fvMesh& mesh)
574 {
575  // Remove any zero-sized ones. Assumes
576  // - processor patches are already only there if needed
577  // - all other patches are available on all processors
578  // - but coupled ones might still be needed, even if zero-size
579  // (e.g. processorCyclic)
580  // See also logic in createPatch.
581  const polyBoundaryMesh& pbm = mesh.boundaryMesh();
582 
583  labelList oldToNew(pbm.size(), -1);
584  label newPatchi = 0;
585  forAll(pbm, patchi)
586  {
587  const polyPatch& pp = pbm[patchi];
588 
589  if (!isA<processorPolyPatch>(pp))
590  {
591  if
592  (
593  isA<coupledPolyPatch>(pp)
594  || returnReduce(pp.size(), sumOp<label>())
595  )
596  {
597  // Coupled (and unknown size) or uncoupled and used
598  oldToNew[patchi] = newPatchi++;
599  }
600  }
601  }
602 
603  forAll(pbm, patchi)
604  {
605  const polyPatch& pp = pbm[patchi];
606 
607  if (isA<processorPolyPatch>(pp))
608  {
609  oldToNew[patchi] = newPatchi++;
610  }
611  }
612 
613 
614  const label nKeepPatches = newPatchi;
615 
616  // Shuffle unused ones to end
617  if (nKeepPatches != pbm.size())
618  {
619  Info<< endl
620  << "Removing zero-sized patches:" << endl << incrIndent;
621 
622  forAll(oldToNew, patchi)
623  {
624  if (oldToNew[patchi] == -1)
625  {
626  Info<< indent << pbm[patchi].name()
627  << " type " << pbm[patchi].type()
628  << " at position " << patchi << endl;
629  oldToNew[patchi] = newPatchi++;
630  }
631  }
632  Info<< decrIndent;
633 
634  fvMeshTools::reorderPatches(mesh, oldToNew, nKeepPatches, true);
635  Info<< endl;
636  }
637 }
638 
639 
640 // Write mesh and additional information
641 void writeMesh
642 (
643  const string& msg,
644  const meshRefinement& meshRefiner,
645  const meshRefinement::debugType debugLevel,
646  const meshRefinement::writeType writeLevel
647 )
648 {
649  const fvMesh& mesh = meshRefiner.mesh();
650 
651  meshRefiner.printMeshInfo(debugLevel, msg);
652  Info<< "Writing mesh to time " << meshRefiner.timeName() << endl;
653 
654  meshRefiner.write
655  (
656  debugLevel,
658  mesh.time().path()/meshRefiner.timeName()
659  );
660  Info<< "Wrote mesh in = "
661  << mesh.time().cpuTimeIncrement() << " s." << endl;
662 }
663 
664 
665 int main(int argc, char *argv[])
666 {
667  #include "addOverwriteOption.H"
669  (
670  "checkGeometry",
671  "check all surface geometry for quality"
672  );
674  (
675  "surfaceSimplify",
676  "boundBox",
677  "simplify the surface using snappyHexMesh starting from a boundBox"
678  );
680  (
681  "patches",
682  "(patch0 .. patchN)",
683  "only triangulate selected patches (wildcards supported)"
684  );
686  (
687  "outFile",
688  "file",
689  "name of the file to save the simplified surface to"
690  );
691  #include "addDictOption.H"
692 
693  #include "setRootCase.H"
694  #include "createTime.H"
695  runTime.functionObjects().off();
696 
697  const bool overwrite = args.optionFound("overwrite");
698  const bool checkGeometry = args.optionFound("checkGeometry");
699  const bool surfaceSimplify = args.optionFound("surfaceSimplify");
700 
702 
703  {
704  Foam::Info
705  << "Create mesh for time = "
706  << runTime.timeName() << Foam::nl << Foam::endl;
707 
708  meshPtr.set
709  (
710  new fvMesh
711  (
713  (
715  runTime.timeName(),
716  runTime,
718  )
719  )
720  );
721  }
722 
723  fvMesh& mesh = meshPtr();
724 
725  Info<< "Read mesh in = "
726  << runTime.cpuTimeIncrement() << " s" << endl;
727 
728  // Check patches and faceZones are synchronised
729  mesh.boundaryMesh().checkParallelSync(true);
731 
732 
733  // Read meshing dictionary
734  const word dictName("snappyHexMeshDict");
735  #include "setSystemMeshDictionaryIO.H"
736  const IOdictionary meshDict(dictIO);
737 
738 
739  // all surface geometry
740  const dictionary& geometryDict = meshDict.subDict("geometry");
741 
742  // refinement parameters
743  const dictionary& refineDict = meshDict.subDict("castellatedMeshControls");
744 
745  // mesh motion and mesh quality parameters
746  const dictionary& motionDict = meshDict.subDict("meshQualityControls");
747 
748  // snap-to-surface parameters
749  const dictionary& snapDict = meshDict.subDict("snapControls");
750 
751  // layer addition parameters
752  const dictionary& layerDict = meshDict.subDict("addLayersControls");
753 
754  // absolute merge distance
755  const scalar mergeDist = getMergeDistance
756  (
757  mesh,
758  meshDict.lookup<scalar>("mergeTolerance")
759  );
760 
761  const Switch keepPatches(meshDict.lookupOrDefault("keepPatches", false));
762 
763 
764 
765  // Read decomposePar dictionary
766  dictionary decomposeDict;
767  {
768  if (Pstream::parRun())
769  {
770  decomposeDict = IOdictionary
771  (
772  IOobject
773  (
774  "decomposeParDict",
775  runTime.system(),
776  mesh,
779  )
780  );
781  }
782  else
783  {
784  decomposeDict.add("method", "none");
785  decomposeDict.add("numberOfSubdomains", 1);
786  }
787  }
788 
789 
790  // Debug
791  // ~~~~~
792 
793  // Set debug level
795  (
796  meshDict.lookupOrDefault<label>
797  (
798  "debug",
799  0
800  )
801  );
802  {
803  wordList flags;
804  if (meshDict.readIfPresent("debugFlags", flags))
805  {
806  debugLevel = meshRefinement::debugType
807  (
809  (
811  flags
812  )
813  );
814  }
815  }
816  if (debugLevel > 0)
817  {
818  meshRefinement::debug = debugLevel;
819  snappyRefineDriver::debug = debugLevel;
820  snappySnapDriver::debug = debugLevel;
821  snappyLayerDriver::debug = debugLevel;
822  }
823 
824  // Set file writing level
825  {
826  wordList flags;
827  if (meshDict.readIfPresent("writeFlags", flags))
828  {
830  (
832  (
834  (
836  flags
837  )
838  )
839  );
840  }
841  }
842 
843  // Set output level
844  {
845  wordList flags;
846  if (meshDict.readIfPresent("outputFlags", flags))
847  {
849  (
851  (
853  (
855  flags
856  )
857  )
858  );
859  }
860  }
861 
862 
863  // Read geometry
864  // ~~~~~~~~~~~~~
865 
866  searchableSurfaces allGeometry
867  (
868  IOobject
869  (
870  "abc", // dummy name
871  mesh.time().constant(), // instance
872  // mesh.time().findInstance("triSurface", word::null),// instance
873  "triSurface", // local
874  mesh.time(), // registry
877  ),
878  geometryDict,
879  meshDict.lookupOrDefault("singleRegionName", true)
880  );
881 
882 
883  // Read refinement surfaces
884  // ~~~~~~~~~~~~~~~~~~~~~~~~
885 
886  autoPtr<refinementSurfaces> surfacesPtr;
887 
888  Info<< "Reading refinement surfaces." << endl;
889 
890  if (surfaceSimplify)
891  {
892  IOdictionary foamyHexMeshDict
893  (
894  IOobject
895  (
896  "foamyHexMeshDict",
897  runTime.system(),
898  runTime,
901  )
902  );
903 
904  const dictionary& conformationDict =
905  foamyHexMeshDict.subDict("surfaceConformation").subDict
906  (
907  "geometryToConformTo"
908  );
909 
910  const dictionary& motionDict =
911  foamyHexMeshDict.subDict("motionControl");
912 
913  const dictionary& shapeControlDict =
914  motionDict.subDict("shapeControlFunctions");
915 
916  // Calculate current ratio of hex cells v.s. wanted cell size
917  const scalar defaultCellSize =
918  motionDict.lookup<scalar>("defaultCellSize");
919 
920  const scalar initialCellSize = ::pow(meshPtr().V()[0], 1.0/3.0);
921 
922  // Info<< "Wanted cell size = " << defaultCellSize << endl;
923  // Info<< "Current cell size = " << initialCellSize << endl;
924  // Info<< "Fraction = " << initialCellSize/defaultCellSize
925  // << endl;
926 
927  surfacesPtr =
928  createRefinementSurfaces
929  (
930  allGeometry,
931  conformationDict,
932  shapeControlDict,
933  refineDict.lookupOrDefault("gapLevelIncrement", 0),
934  initialCellSize/defaultCellSize
935  );
936  }
937  else
938  {
939  surfacesPtr.set
940  (
942  (
943  allGeometry,
944  refineDict.subDict("refinementSurfaces"),
945  refineDict.lookupOrDefault("gapLevelIncrement", 0)
946  )
947  );
948 
949  Info<< "Read refinement surfaces in = "
950  << mesh.time().cpuTimeIncrement() << " s" << nl << endl;
951  }
952 
953  refinementSurfaces& surfaces = surfacesPtr();
954 
955 
956  // Checking only?
957 
958  if (checkGeometry)
959  {
960  // Extract patchInfo
961  List<wordList> patchTypes(allGeometry.size());
962 
963  const PtrList<dictionary>& patchInfo = surfaces.patchInfo();
964  const labelList& surfaceGeometry = surfaces.surfaces();
965  forAll(surfaceGeometry, surfI)
966  {
967  label geomI = surfaceGeometry[surfI];
968  const wordList& regNames = allGeometry.regionNames()[geomI];
969 
970  patchTypes[geomI].setSize(regNames.size());
971  forAll(regNames, regionI)
972  {
973  label globalRegionI = surfaces.globalRegion(surfI, regionI);
974 
975  if (patchInfo.set(globalRegionI))
976  {
977  patchTypes[geomI][regionI] =
978  word(patchInfo[globalRegionI].lookup("type"));
979  }
980  else
981  {
982  patchTypes[geomI][regionI] = wallPolyPatch::typeName;
983  }
984  }
985  }
986 
987  // Write some stats
988  allGeometry.writeStats(patchTypes, Info);
989  // Check topology
990  allGeometry.checkTopology(true);
991  // Check geometry
992  allGeometry.checkGeometry
993  (
994  100.0, // max size ratio
995  1e-9, // intersection tolerance
997  0.01, // min triangle quality
998  true
999  );
1000 
1001  return 0;
1002  }
1003 
1004 
1005 
1006  // Read refinement shells
1007  // ~~~~~~~~~~~~~~~~~~~~~~
1008 
1009  Info<< "Reading refinement shells." << endl;
1010  shellSurfaces shells
1011  (
1012  allGeometry,
1013  refineDict.subDict("refinementRegions")
1014  );
1015  Info<< "Read refinement shells in = "
1016  << mesh.time().cpuTimeIncrement() << " s" << nl << endl;
1017 
1018 
1019  Info<< "Setting refinement level of surface to be consistent"
1020  << " with shells." << endl;
1021  surfaces.setMinLevelFields(shells);
1022  Info<< "Checked shell refinement in = "
1023  << mesh.time().cpuTimeIncrement() << " s" << nl << endl;
1024 
1025 
1026  // Read feature meshes
1027  // ~~~~~~~~~~~~~~~~~~~
1028 
1029  Info<< "Reading features." << endl;
1030  refinementFeatures features
1031  (
1032  mesh,
1033  refineDict.lookup("features")
1034  );
1035  Info<< "Read features in = "
1036  << mesh.time().cpuTimeIncrement() << " s" << nl << endl;
1037 
1038 
1039 
1040  // Refinement engine
1041  // ~~~~~~~~~~~~~~~~~
1042 
1043  Info<< nl
1044  << "Determining initial surface intersections" << nl
1045  << "-----------------------------------------" << nl
1046  << endl;
1047 
1048  // Main refinement engine
1049  meshRefinement meshRefiner
1050  (
1051  mesh,
1052  mergeDist, // tolerance used in sorting coordinates
1053  overwrite, // overwrite mesh files?
1054  surfaces, // for surface intersection refinement
1055  features, // for feature edges/point based refinement
1056  shells // for volume (inside/outside) refinement
1057  );
1058  Info<< "Calculated surface intersections in = "
1059  << mesh.time().cpuTimeIncrement() << " s" << nl << endl;
1060 
1061  // Some stats
1062  meshRefiner.printMeshInfo(debugLevel, "Initial mesh");
1063 
1064  meshRefiner.write
1065  (
1066  meshRefinement::debugType(debugLevel&meshRefinement::OBJINTERSECTIONS),
1068  mesh.time().path()/meshRefiner.timeName()
1069  );
1070 
1071 
1072  // Add all the surface regions as patches
1073  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
1074 
1075  //- Global surface region to patch (non faceZone surface) or patches
1076  // (faceZone surfaces)
1077  labelList globalToMasterPatch;
1078  labelList globalToSlavePatch;
1079  {
1080  Info<< nl
1081  << "Adding patches for surface regions" << nl
1082  << "----------------------------------" << nl
1083  << endl;
1084 
1085  // From global region number to mesh patch.
1086  globalToMasterPatch.setSize(surfaces.nRegions(), -1);
1087  globalToSlavePatch.setSize(surfaces.nRegions(), -1);
1088 
1089  Info<< setf(ios_base::left)
1090  << setw(6) << "Patch"
1091  << setw(20) << "Type"
1092  << setw(30) << "Region" << nl
1093  << setw(6) << "-----"
1094  << setw(20) << "----"
1095  << setw(30) << "------" << endl;
1096 
1097  const labelList& surfaceGeometry = surfaces.surfaces();
1098  const PtrList<dictionary>& surfacePatchInfo = surfaces.patchInfo();
1099 
1100  forAll(surfaceGeometry, surfI)
1101  {
1102  label geomI = surfaceGeometry[surfI];
1103 
1104  const wordList& regNames = allGeometry.regionNames()[geomI];
1105 
1106  Info<< surfaces.names()[surfI] << ':' << nl << nl;
1107 
1108  if (surfaces.surfZones()[surfI].faceZoneName().empty())
1109  {
1110  // 'Normal' surface
1111  forAll(regNames, i)
1112  {
1113  label globalRegionI = surfaces.globalRegion(surfI, i);
1114 
1115  label patchi;
1116 
1117  if (surfacePatchInfo.set(globalRegionI))
1118  {
1119  patchi = meshRefiner.addMeshedPatch
1120  (
1121  regNames[i],
1122  surfacePatchInfo[globalRegionI]
1123  );
1124  }
1125  else
1126  {
1127  dictionary patchInfo;
1128  patchInfo.set("type", wallPolyPatch::typeName);
1129 
1130  patchi = meshRefiner.addMeshedPatch
1131  (
1132  regNames[i],
1133  patchInfo
1134  );
1135  }
1136 
1137  Info<< setf(ios_base::left)
1138  << setw(6) << patchi
1139  << setw(20) << mesh.boundaryMesh()[patchi].type()
1140  << setw(30) << regNames[i] << nl;
1141 
1142  globalToMasterPatch[globalRegionI] = patchi;
1143  globalToSlavePatch[globalRegionI] = patchi;
1144  }
1145  }
1146  else
1147  {
1148  // Zoned surface
1149  forAll(regNames, i)
1150  {
1151  label globalRegionI = surfaces.globalRegion(surfI, i);
1152 
1153  // Add master side patch
1154  {
1155  label patchi;
1156 
1157  if (surfacePatchInfo.set(globalRegionI))
1158  {
1159  patchi = meshRefiner.addMeshedPatch
1160  (
1161  regNames[i],
1162  surfacePatchInfo[globalRegionI]
1163  );
1164  }
1165  else
1166  {
1167  dictionary patchInfo;
1168  patchInfo.set("type", wallPolyPatch::typeName);
1169 
1170  patchi = meshRefiner.addMeshedPatch
1171  (
1172  regNames[i],
1173  patchInfo
1174  );
1175  }
1176 
1177  Info<< setf(ios_base::left)
1178  << setw(6) << patchi
1179  << setw(20) << mesh.boundaryMesh()[patchi].type()
1180  << setw(30) << regNames[i] << nl;
1181 
1182  globalToMasterPatch[globalRegionI] = patchi;
1183  }
1184  // Add slave side patch
1185  {
1186  const word slaveName = regNames[i] + "_slave";
1187  label patchi;
1188 
1189  if (surfacePatchInfo.set(globalRegionI))
1190  {
1191  patchi = meshRefiner.addMeshedPatch
1192  (
1193  slaveName,
1194  surfacePatchInfo[globalRegionI]
1195  );
1196  }
1197  else
1198  {
1199  dictionary patchInfo;
1200  patchInfo.set("type", wallPolyPatch::typeName);
1201 
1202  patchi = meshRefiner.addMeshedPatch
1203  (
1204  slaveName,
1205  patchInfo
1206  );
1207  }
1208 
1209  Info<< setf(ios_base::left)
1210  << setw(6) << patchi
1211  << setw(20) << mesh.boundaryMesh()[patchi].type()
1212  << setw(30) << slaveName << nl;
1213 
1214  globalToSlavePatch[globalRegionI] = patchi;
1215  }
1216  }
1217  }
1218 
1219  Info<< nl;
1220  }
1221  Info<< "Added patches in = "
1222  << mesh.time().cpuTimeIncrement() << " s" << nl << endl;
1223  }
1224 
1225 
1226  // Parallel
1227  // ~~~~~~~~
1228 
1229  // Decomposition
1230  autoPtr<decompositionMethod> decomposerPtr
1231  (
1233  (
1234  decomposeDict
1235  )
1236  );
1237  decompositionMethod& decomposer = decomposerPtr();
1238 
1239  if (Pstream::parRun() && !decomposer.parallelAware())
1240  {
1242  << "You have selected decomposition method "
1243  << decomposer.typeName
1244  << " which is not parallel aware." << endl
1245  << "Please select one that is (hierarchical, ptscotch)"
1246  << exit(FatalError);
1247  }
1248 
1249  // Mesh distribution engine (uses tolerance to reconstruct meshes)
1250  fvMeshDistribute distributor(mesh, mergeDist);
1251 
1252 
1253 
1254 
1255 
1256  // Now do the real work -refinement -snapping -layers
1257  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
1258 
1259  const Switch wantRefine(meshDict.lookup("castellatedMesh"));
1260  const Switch wantSnap(meshDict.lookup("snap"));
1261  const Switch wantLayers(meshDict.lookup("addLayers"));
1262 
1263  // Refinement parameters
1264  const refinementParameters refineParams(refineDict);
1265 
1266  // Snap parameters
1267  const snapParameters snapParams(snapDict);
1268 
1269  // Layer addition parameters
1270  const layerParameters layerParams(layerDict, mesh.boundaryMesh());
1271 
1272 
1273  if (wantRefine)
1274  {
1275  cpuTime timer;
1276 
1277  snappyRefineDriver refineDriver
1278  (
1279  meshRefiner,
1280  decomposer,
1281  distributor,
1282  globalToMasterPatch,
1283  globalToSlavePatch
1284  );
1285 
1286 
1287  if (!overwrite && !debugLevel)
1288  {
1289  const_cast<Time&>(mesh.time())++;
1290  }
1291 
1292 
1293  refineDriver.doRefine
1294  (
1295  refineDict,
1296  refineParams,
1297  snapParams,
1298  refineParams.handleSnapProblems(),
1299  motionDict
1300  );
1301 
1302 
1303  if (!keepPatches && !wantSnap && !wantLayers)
1304  {
1305  removeZeroSizedPatches(mesh);
1306  }
1307 
1308  writeMesh
1309  (
1310  "Refined mesh",
1311  meshRefiner,
1312  debugLevel,
1314  );
1315 
1316  Info<< "Mesh refined in = "
1317  << timer.cpuTimeIncrement() << " s." << endl;
1318  }
1319 
1320  if (wantSnap)
1321  {
1322  cpuTime timer;
1323 
1324  snappySnapDriver snapDriver
1325  (
1326  meshRefiner,
1327  globalToMasterPatch,
1328  globalToSlavePatch
1329  );
1330 
1331  if (!overwrite && !debugLevel)
1332  {
1333  const_cast<Time&>(mesh.time())++;
1334  }
1335 
1336  // Use the resolveFeatureAngle from the refinement parameters
1337  scalar curvature = refineParams.curvature();
1338  scalar planarAngle = refineParams.planarAngle();
1339 
1340  snapDriver.doSnap
1341  (
1342  snapDict,
1343  motionDict,
1344  curvature,
1345  planarAngle,
1346  snapParams
1347  );
1348 
1349  if (!keepPatches && !wantLayers)
1350  {
1351  removeZeroSizedPatches(mesh);
1352  }
1353 
1354  writeMesh
1355  (
1356  "Snapped mesh",
1357  meshRefiner,
1358  debugLevel,
1360  );
1361 
1362  Info<< "Mesh snapped in = "
1363  << timer.cpuTimeIncrement() << " s." << endl;
1364  }
1365 
1366  if (wantLayers)
1367  {
1368  cpuTime timer;
1369 
1370  snappyLayerDriver layerDriver
1371  (
1372  meshRefiner,
1373  globalToMasterPatch,
1374  globalToSlavePatch
1375  );
1376 
1377  // Use the maxLocalCells from the refinement parameters
1378  bool preBalance = returnReduce
1379  (
1380  (mesh.nCells() >= refineParams.maxLocalCells()),
1381  orOp<bool>()
1382  );
1383 
1384 
1385  if (!overwrite && !debugLevel)
1386  {
1387  const_cast<Time&>(mesh.time())++;
1388  }
1389 
1390  layerDriver.doLayers
1391  (
1392  layerDict,
1393  motionDict,
1394  layerParams,
1395  preBalance,
1396  decomposer,
1397  distributor
1398  );
1399 
1400  if (!keepPatches)
1401  {
1402  removeZeroSizedPatches(mesh);
1403  }
1404 
1405  writeMesh
1406  (
1407  "Layer mesh",
1408  meshRefiner,
1409  debugLevel,
1411  );
1412 
1413  Info<< "Layers added in = "
1414  << timer.cpuTimeIncrement() << " s." << endl;
1415  }
1416 
1417 
1418  {
1419  // Check final mesh
1420  Info<< "Checking final mesh ..." << endl;
1421  faceSet wrongFaces(mesh, "wrongFaces", mesh.nFaces()/100);
1422  motionSmoother::checkMesh(false, mesh, motionDict, wrongFaces);
1423  const label nErrors = returnReduce
1424  (
1425  wrongFaces.size(),
1426  sumOp<label>()
1427  );
1428 
1429  if (nErrors > 0)
1430  {
1431  Info<< "Finished meshing with " << nErrors << " illegal faces"
1432  << " (concave, zero area or negative cell pyramid volume)"
1433  << endl;
1434  wrongFaces.write();
1435  }
1436  else
1437  {
1438  Info<< "Finished meshing without any errors" << endl;
1439  }
1440  }
1441 
1442 
1443  if (surfaceSimplify)
1444  {
1445  const polyBoundaryMesh& bMesh = mesh.boundaryMesh();
1446 
1447  labelHashSet includePatches(bMesh.size());
1448 
1449  if (args.optionFound("patches"))
1450  {
1451  includePatches = bMesh.patchSet
1452  (
1453  wordReList(args.optionLookup("patches")())
1454  );
1455  }
1456  else
1457  {
1458  forAll(bMesh, patchi)
1459  {
1460  const polyPatch& patch = bMesh[patchi];
1461 
1462  if (!isA<processorPolyPatch>(patch))
1463  {
1464  includePatches.insert(patchi);
1465  }
1466  }
1467  }
1468 
1469  fileName outFileName
1470  (
1472  (
1473  "outFile",
1474  "constant/triSurface/simplifiedSurface.stl"
1475  )
1476  );
1477 
1478  extractSurface
1479  (
1480  mesh,
1481  runTime,
1482  includePatches,
1483  outFileName
1484  );
1485 
1486  pointIOField cellCentres
1487  (
1488  IOobject
1489  (
1490  "internalCellCentres",
1491  runTime.timeName(),
1492  mesh,
1495  ),
1496  mesh.cellCentres()
1497  );
1498 
1499  cellCentres.write();
1500  }
1501 
1502 
1503  Info<< "Finished meshing in = "
1504  << runTime.elapsedCpuTime() << " s." << endl;
1505 
1506  Info<< "End\n" << endl;
1507 
1508  return 0;
1509 }
1510 
1511 
1512 // ************************************************************************* //
static bool checkMesh(const bool report, const polyMesh &mesh, const dictionary &dict, labelHashSet &wrongFaces)
Check mesh with mesh settings in dict. Collects incorrect faces.
A surface geometry mesh, in which the surface zone information is conveyed by the &#39;zoneId&#39; associated...
Definition: MeshedSurface.H:74
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:434
void printMeshInfo(const bool, const string &) const
Print some mesh stats.
A HashTable with keys but without contents.
Definition: HashSet.H:59
static void reorderPatches(fvMesh &, const labelList &oldToNew, const label nPatches, const bool validBoundary)
Reorder and remove trailing patches. If validBoundary call is parallel.
Definition: fvMeshTools.C:140
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:667
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
bool processorCase() const
Return true if this is a processor case.
Definition: TimePaths.H:90
fileName path() const
Return path.
Definition: Time.H:265
Simple container to keep together layer specific information.
virtual Ostream & write(const char)=0
Write character.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
const word & name() const
Return name.
const word & name() const
Return name.
Definition: IOobject.H:303
A class for handling file names.
Definition: fileName.H:79
const labelList & surfaces() const
Ostream & indent(Ostream &os)
Indent stream.
Definition: Ostream.H:221
A list of face labels.
Definition: faceSet.H:48
virtual bool parallelAware() const =0
Is method parallel aware (i.e. does it synchronize domains across.
const entry * lookupEntryPtr(const word &, bool recursive, bool patternMatch) const
Find and return an entry data stream pointer if present.
Definition: dictionary.C:706
bool set(const label) const
Is element set.
Definition: PtrListI.H:65
An identifier for a surface zone on a meshed surface.
void inplaceRenumber(const labelUList &oldToNew, ListType &)
Inplace renumber the values of a list.
const keyType & keyword() const
Return keyword.
Definition: entry.H:123
label checkTopology(const bool report) const
All topological checks. Return number of failed checks.
dimensionedScalar log(const dimensionedScalar &ds)
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
void off()
Switch the function objects off.
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:158
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
bool clean()
Cleanup file name.
Definition: fileName.C:81
label nFaces() const
static fvMesh * meshPtr
Definition: globalFoam.H:52
static unsigned int defaultPrecision()
Return the default precision.
Definition: IOstream.H:461
static const NamedEnum< IOoutputType, 1 > IOoutputTypeNames
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
wordList patchTypes(nPatches)
static int myProcNo(const label communicator=0)
Number of this process (starting from masterNo() = 0)
Definition: UPstream.H:429
label nCells() const
engineTime & runTime
static word defaultRegion
Return the default region name.
Definition: polyMesh.H:309
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
bool optionFound(const word &opt) const
Return true if the named option is found.
Definition: argListI.H:114
A simple wrapper around bool so that it can be read as a word: true/false, on/off, yes/no, y/n, t/f, or none/any.
Definition: Switch.H:60
static bool master(const label communicator=0)
Am I the master process.
Definition: UPstream.H:423
static writeType writeLevel()
Get/set write level.
virtual const dictionary & dict() const =0
Return dictionary if this entry is a dictionary.
label findPatchID(const word &patchName) const
Find patch index given a name.
A bounding box defined in terms of the points at its extremities.
Definition: boundBox.H:58
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:622
Base class of (analytical or triangulated) surface. Encapsulates all the search routines. WIP.
autoPtr< globalIndex > mergePoints(labelList &pointToGlobal, labelList &uniquePoints) const
Helper for merging (collocated!) mesh point data.
const labelList & minLevel() const
From global region number to refinement level.
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:239
Simple container to keep together refinement specific information.
bool insert(const Key &key)
Insert a new entry.
Definition: HashSet.H:111
labelHashSet patchSet(const UList< wordRe > &patchNames, const bool warnNotFound=true, const bool usePatchGroups=true) const
Return the set of patch IDs corresponding to the given names.
bool erase(T *p)
Remove the specified element from the list and delete it.
Definition: ILList.C:105
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:68
autoPtr< dictionary > clone() const
Construct and return clone.
Definition: dictionary.C:576
wordList toc() const
Return the table of contents.
Definition: dictionary.C:1018
const fvMesh & mesh() const
Reference to mesh.
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Definition: IOdictionary.H:52
Container for data on surfaces used for surface-driven refinement. Contains all the data about the le...
virtual const pointField & points() const
Return raw points.
Definition: polyMesh.C:1131
static void mapCombineScatter(const List< commsStruct > &comms, Container &Values, const int tag, const label comm)
Scatter data. Reverse of combineGather.
const dictionary & subDict(const word &) const
Find and return a sub-dictionary.
Definition: dictionary.C:934
Encapsulates queries for features.
T optionLookupOrDefault(const word &opt, const T &deflt) const
Read a value from the named option if present.
Definition: argListI.H:243
A list of faces which address into the list of points.
Sends/receives parts of mesh+fvfields to neighbouring processors. Used in load balancing.
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:42
stressControl lookup("compactNormalStress") >> compactNormalStress
T clone(const T &t)
Definition: List.H:54
dynamicFvMesh & mesh
const dictionary & optionalSubDict(const word &) const
Find and return a sub-dictionary if found.
Definition: dictionary.C:1001
A class for handling words, derived from string.
Definition: word.H:59
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Simple container to keep together snap specific information.
double cpuTimeIncrement() const
Return CPU time (in seconds) since last call to cpuTimeIncrement()
Definition: cpuTime.C:74
const labelList & gapLevel() const
From global region number to small gap refinement level.
const List< wordList > & regionNames() const
scalar mag() const
The magnitude of the bounding box span.
Definition: boundBoxI.H:90
const word & constant() const
Return constant name.
Definition: TimePaths.H:124
static void addOption(const word &opt, const string &param="", const string &usage="")
Add to an option to validOptions with usage information.
Definition: argList.C:128
static const NamedEnum< IOdebugType, 5 > IOdebugTypeNames
static void checkCoupledFaceZones(const polyMesh &)
Helper function: check that face zones are synced.
Container for searchableSurfaces.
Smanip< ios_base::fmtflags > setf(const ios_base::fmtflags flags)
Definition: IOmanip.H:164
const globalMeshData & globalData() const
Return parallel info.
Definition: polyMesh.C:1394
const word dictName("particleTrackDict")
Abstract base class for decomposition.
Encapsulates queries for volume refinement (&#39;refine all cells within shell&#39;).
Definition: shellSurfaces.H:53
static void mapCombineGather(const List< commsStruct > &comms, Container &Values, const CombineOp &cop, const int tag, const label comm)
virtual const faceList & faces() const
Return raw faces.
Definition: polyMesh.C:1156
const PtrList< surfaceZonesInfo > & surfZones() const
const vectorField & cellCentres() const
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
word timeName() const
Replacement for Time::timeName() : return oldInstance (if.
const word & system() const
Return system name.
Definition: TimePaths.H:114
static autoPtr< decompositionMethod > New(const dictionary &decompositionDict)
Return a reference to the selected decomposition method.
Foam::polyBoundaryMesh.
Istream and Ostream manipulators taking arguments.
bool checkParallelSync(const bool report=false) const
Check whether all procs have all patches and in same order. Return.
static const char nl
Definition: Ostream.H:260
IOstream::streamFormat writeFormat() const
Default write format.
Definition: Time.H:287
const Time & time() const
Return time.
static int readFlags(const Enum &namedEnum, const wordList &)
Helper: convert wordList into bit pattern using provided.
label addMeshedPatch(const word &name, const dictionary &)
Add patch originating from meshing. Update meshedPatches_.
void set(T *)
Set pointer to that given.
Definition: autoPtrI.H:99
Ostream & decrIndent(Ostream &os)
Decrement the indent level.
Definition: Ostream.H:235
labelHashSet includePatches
All to do with adding layers.
static outputType outputLevel()
Get/set output level.
IOobject dictIO(dictName, runTime.constant(), mesh, IOobject::MUST_READ_IF_MODIFIED, IOobject::NO_WRITE)
Helper class which maintains intersections of (changing) mesh with (static) surfaces.
label size() const
Return the number of elements in the UPtrList.
Definition: UPtrListI.H:29
label checkGeometry(const scalar maxRatio, const scalar tolerance, const autoPtr< writer< scalar >> &setWriter, const scalar minQuality, const bool report) const
All geometric checks. Return number of failed checks.
All to do with snapping to surface.
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
void setSize(const label)
Reset size of List.
Definition: List.C:281
const boundBox & bounds() const
Return mesh bounding box.
Definition: polyMesh.H:440
static bool & parRun()
Is this a parallel run?
Definition: UPstream.H:399
const wordList & names() const
Names of surfaces.
static label nProcs(const label communicator=0)
Number of processes in parallel run.
Definition: UPstream.H:411
label patchi
double elapsedCpuTime() const
Return CPU time (in seconds) from the start.
Definition: cpuTime.C:67
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
const functionObjectList & functionObjects() const
Return the list of function objects.
Definition: Time.H:409
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:70
const wordList & names() const
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
A List with indirect addressing.
Definition: fvMatrix.H:106
label index() const
Return the index of this patch in the boundaryMesh.
label start() const
Return start label of this patch in the polyMesh face list.
Definition: polyPatch.H:303
label globalRegion(const label surfI, const label regionI) const
From surface and region on surface to global region.
void set(entry *)
Assign a new entry, overwrite any existing entry.
Definition: dictionary.C:1194
messageStream Info
const PtrList< dictionary > & patchInfo() const
From global region number to patch type.
tmp< pointField > allPoints(const Triangulation &t)
Extract all points in vertex-index order.
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:52
const doubleScalar e
Elementary charge.
Definition: doubleScalar.H:105
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
static void addBoolOption(const word &opt, const string &usage="")
Add to a bool option to validOptions with usage information.
Definition: argList.C:118
Omanip< int > setw(const int i)
Definition: IOmanip.H:199
Starts timing CPU usage and return elapsed time from start.
Definition: cpuTime.H:52
List< wordRe > wordReList
A List of wordRe (word or regular expression)
Definition: wordReList.H:50
bool write() const
Write mesh and all data.
virtual bool write(const bool write=true) const
Write using setting from DB.
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
void setMinLevelFields(const shellSurfaces &shells)
Calculate minLevelFields.
Foam::argList args(argc, argv)
void writeStats(const List< wordList > &, Ostream &) const
Write some stats.
Ostream & incrIndent(Ostream &os)
Increment the indent level.
Definition: Ostream.H:228
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
static void gatherList(const List< commsStruct > &comms, List< T > &Values, const int tag, const label comm)
Gather data but keep individual values separate.
A primitive field of type <T> with automated input and output.
Definition: IOField.H:50
Namespace for OpenFOAM.
A keyword and a list of tokens is an &#39;entry&#39;.
Definition: entry.H:65
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:812
static const NamedEnum< IOwriteType, 5 > IOwriteTypeNames
const labelList & maxLevel() const
From global region number to refinement level.
IStringStream optionLookup(const word &opt) const
Return an IStringStream from the named option.
Definition: argListI.H:120