All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
blockMesh.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-2018 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  blockMesh
26 
27 Description
28  A multi-block mesh generator.
29 
30  Uses the block mesh description found in
31  - \c system/blockMeshDict
32  - \c system/<region>/blockMeshDict
33  - \c constant/polyMesh/blockMeshDict
34  - \c constant/<region>/polyMesh/blockMeshDict
35 
36 Usage
37  \b blockMesh [OPTION]
38 
39  Options:
40  - \par -blockTopology
41  Write the topology as a set of edges in OBJ format.
42 
43  - \par -region <name>
44  Specify an alternative mesh region.
45 
46  - \par -dict <filename>
47  Specify alternative dictionary for the block mesh description.
48 
49 \*---------------------------------------------------------------------------*/
50 
51 #include "Time.H"
52 #include "IOdictionary.H"
53 #include "IOPtrList.H"
54 
55 #include "blockMesh.H"
56 #include "attachPolyTopoChanger.H"
57 #include "polyTopoChange.H"
58 #include "emptyPolyPatch.H"
59 #include "cyclicPolyPatch.H"
60 
61 #include "argList.H"
62 #include "OSspecific.H"
63 #include "OFstream.H"
64 
65 #include "Pair.H"
66 #include "slidingInterface.H"
67 
68 using namespace Foam;
69 
70 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
71 
72 int main(int argc, char *argv[])
73 {
76  (
77  "blockTopology",
78  "write block edges and centres as .obj files"
79  );
81  (
82  "noClean",
83  "keep the existing files in the polyMesh"
84  );
86  (
87  "dict",
88  "file",
89  "specify alternative dictionary for the blockMesh description"
90  );
91 
93  (
94  "Block description\n"
95  "\n"
96  " For a given block, the correspondence between the ordering of\n"
97  " vertex labels and face labels is shown below.\n"
98  " For vertex numbering in the sequence 0 to 7 (block, centre):\n"
99  " faces 0 (f0) and 1 are left and right, respectively;\n"
100  " faces 2 and 3 are front and back; \n"
101  " and faces 4 and 5 are bottom and top::\n"
102  "\n"
103  " 7 ---- 6\n"
104  " f5 |\\ |\\ f3\n"
105  " | | 4 ---- 5 \\\n"
106  " | 3 |--- 2 | \\\n"
107  " | \\| \\| f2\n"
108  " f4 0 ---- 1\n"
109  "\n"
110  " Z f0 ----- f1\n"
111  " | Y\n"
112  " | /\n"
113  " O --- X\n"
114  );
115 
116  #include "addRegionOption.H"
117  #include "setRootCase.H"
118  #include "createTime.H"
119 
120  const word dictName("blockMeshDict");
121 
123  word regionPath;
124 
125  // Check if the region is specified otherwise mesh the default region
126  if (args.optionReadIfPresent("region", regionName, polyMesh::defaultRegion))
127  {
128  Info<< nl << "Generating mesh for region " << regionName << endl;
129  regionPath = regionName;
130  }
131 
132  // Search for the appropriate blockMesh dictionary....
133 
135 
136  // Check if the dictionary is specified on the command-line
137  if (args.optionFound("dict"))
138  {
139  dictPath = args["dict"];
140 
141  dictPath =
142  (
143  isDir(dictPath)
144  ? dictPath/dictName
145  : dictPath
146  );
147  }
148  // Check if dictionary is present in the constant directory
149  else if
150  (
151  exists
152  (
154  /regionPath/polyMesh::meshSubDir/dictName
155  )
156  )
157  {
158  dictPath =
159  runTime.constant()
160  /regionPath/polyMesh::meshSubDir/dictName;
161  }
162  // Otherwise assume the dictionary is present in the system directory
163  else
164  {
165  dictPath = runTime.system()/regionPath/dictName;
166  }
167 
168  if (!args.optionFound("noClean"))
169  {
170  fileName polyMeshPath
171  (
173  );
174 
175  if (exists(polyMeshPath))
176  {
177  if (exists(polyMeshPath/dictName))
178  {
179  Info<< "Not deleting polyMesh directory " << nl
180  << " " << polyMeshPath << nl
181  << " because it contains " << dictName << endl;
182  }
183  else
184  {
185  Info<< "Deleting polyMesh directory" << nl
186  << " " << polyMeshPath << endl;
187  rmDir(polyMeshPath);
188  }
189  }
190  }
191 
192  IOobject meshDictIO
193  (
194  dictPath,
195  runTime,
198  false
199  );
200 
201  if (!meshDictIO.typeHeaderOk<IOdictionary>(true))
202  {
204  << meshDictIO.objectPath()
205  << nl
206  << exit(FatalError);
207  }
208 
209  Info<< "Creating block mesh from\n "
210  << meshDictIO.objectPath() << endl;
211 
212  IOdictionary meshDict(meshDictIO);
213  blockMesh blocks(meshDict, regionName);
214 
215 
216  if (args.optionFound("blockTopology"))
217  {
218  // Write mesh as edges.
219  {
220  fileName objMeshFile("blockTopology.obj");
221 
222  OFstream str(runTime.path()/objMeshFile);
223 
224  Info<< nl << "Dumping block structure as Lightwave obj format"
225  << " to " << objMeshFile << endl;
226 
227  blocks.writeTopology(str);
228  }
229 
230  // Write centres of blocks
231  {
232  fileName objCcFile("blockCentres.obj");
233 
234  OFstream str(runTime.path()/objCcFile);
235 
236  Info<< nl << "Dumping block centres as Lightwave obj format"
237  << " to " << objCcFile << endl;
238 
239  const polyMesh& topo = blocks.topology();
240 
241  const pointField& cellCentres = topo.cellCentres();
242 
243  forAll(cellCentres, celli)
244  {
245  // point cc = b.blockShape().centre(b.points());
246  const point& cc = cellCentres[celli];
247 
248  str << "v " << cc.x() << ' ' << cc.y() << ' ' << cc.z() << nl;
249  }
250  }
251 
252  Info<< nl << "end" << endl;
253 
254  return 0;
255  }
256 
257 
258  Info<< nl << "Creating polyMesh from blockMesh" << endl;
259 
260  word defaultFacesName = "defaultFaces";
261  word defaultFacesType = emptyPolyPatch::typeName;
262  polyMesh mesh
263  (
264  IOobject
265  (
266  regionName,
267  runTime.constant(),
268  runTime
269  ),
270  xferCopy(blocks.points()), // could we re-use space?
271  blocks.cells(),
272  blocks.patches(),
273  blocks.patchNames(),
274  blocks.patchDicts(),
277  );
278 
279 
280  // Read in a list of dictionaries for the merge patch pairs
281  if (meshDict.found("mergePatchPairs"))
282  {
283  List<Pair<word>> mergePatchPairs
284  (
285  meshDict.lookup("mergePatchPairs")
286  );
287 
288  #include "mergePatchPairs.H"
289  }
290  else
291  {
292  Info<< nl << "There are no merge patch pairs edges" << endl;
293  }
294 
295 
296  // Set any cellZones (note: cell labelling unaffected by above
297  // mergePatchPairs)
298 
299  label nZones = blocks.numZonedBlocks();
300 
301  if (nZones > 0)
302  {
303  Info<< nl << "Adding cell zones" << endl;
304 
305  // Map from zoneName to cellZone index
306  HashTable<label> zoneMap(nZones);
307 
308  // Cells per zone.
309  List<DynamicList<label>> zoneCells(nZones);
310 
311  // Running cell counter
312  label celli = 0;
313 
314  // Largest zone so far
315  label freeZoneI = 0;
316 
317  forAll(blocks, blockI)
318  {
319  const block& b = blocks[blockI];
320  const List<FixedList<label, 8>> blockCells = b.cells();
321  const word& zoneName = b.zoneName();
322 
323  if (zoneName.size())
324  {
325  HashTable<label>::const_iterator iter = zoneMap.find(zoneName);
326 
327  label zoneI;
328 
329  if (iter == zoneMap.end())
330  {
331  zoneI = freeZoneI++;
332 
333  Info<< " " << zoneI << '\t' << zoneName << endl;
334 
335  zoneMap.insert(zoneName, zoneI);
336  }
337  else
338  {
339  zoneI = iter();
340  }
341 
342  forAll(blockCells, i)
343  {
344  zoneCells[zoneI].append(celli++);
345  }
346  }
347  else
348  {
349  celli += b.cells().size();
350  }
351  }
352 
353 
354  List<cellZone*> cz(zoneMap.size());
355 
356  forAllConstIter(HashTable<label>, zoneMap, iter)
357  {
358  label zoneI = iter();
359 
360  cz[zoneI] = new cellZone
361  (
362  iter.key(),
363  zoneCells[zoneI].shrink(),
364  zoneI,
365  mesh.cellZones()
366  );
367  }
368 
369  mesh.pointZones().setSize(0);
370  mesh.faceZones().setSize(0);
371  mesh.cellZones().setSize(0);
373  }
374 
375 
376  // Detect any cyclic patches and force re-ordering of the faces
377  {
378  const polyPatchList& patches = mesh.boundaryMesh();
379  bool hasCyclic = false;
380  forAll(patches, patchi)
381  {
382  if (isA<cyclicPolyPatch>(patches[patchi]))
383  {
384  hasCyclic = true;
385  break;
386  }
387  }
388 
389  if (hasCyclic)
390  {
391  Info<< nl << "Detected cyclic patches; ordering boundary faces"
392  << endl;
393  const word oldInstance = mesh.instance();
394  polyTopoChange meshMod(mesh);
395  meshMod.changeMesh(mesh, false);
396  mesh.setInstance(oldInstance);
397  }
398  }
399 
400 
401  // Set the precision of the points data to 10
403 
404  Info<< nl << "Writing polyMesh" << endl;
405  mesh.removeFiles();
406  if (!mesh.write())
407  {
409  << "Failed writing polyMesh."
410  << exit(FatalError);
411  }
412 
413 
414  // Write summary
415  {
416  const polyPatchList& patches = mesh.boundaryMesh();
417 
418  Info<< "----------------" << nl
419  << "Mesh Information" << nl
420  << "----------------" << nl
421  << " " << "boundingBox: " << boundBox(mesh.points()) << nl
422  << " " << "nPoints: " << mesh.nPoints() << nl
423  << " " << "nCells: " << mesh.nCells() << nl
424  << " " << "nFaces: " << mesh.nFaces() << nl
425  << " " << "nInternalFaces: " << mesh.nInternalFaces() << nl;
426 
427  Info<< "----------------" << nl
428  << "Patches" << nl
429  << "----------------" << nl;
430 
431  forAll(patches, patchi)
432  {
433  const polyPatch& p = patches[patchi];
434 
435  Info<< " " << "patch " << patchi
436  << " (start: " << p.start()
437  << " size: " << p.size()
438  << ") name: " << p.name()
439  << nl;
440  }
441  }
442 
443  Info<< "\nEnd\n" << endl;
444 
445  return 0;
446 }
447 
448 
449 // ************************************************************************* //
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:424
word defaultFacesType
Definition: readKivaGrid.H:461
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
fileName path() const
Return path.
Definition: Time.H:266
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.
A class for handling file names.
Definition: fileName.H:69
const faceZoneMesh & faceZones() const
Return face zone mesh.
Definition: polyMesh.H:466
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
static iteratorEnd end()
iteratorEnd set to beyond the end of any HashTable
Definition: HashTable.H:110
error FatalError
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
label nInternalFaces() const
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: HashTable.H:60
virtual bool write(const bool valid=true) const
Write mesh using IO settings from time.
Definition: fvMesh.C:877
label nFaces() const
Foam::word regionName
Output to file stream.
Definition: OFstream.H:82
static word meshSubDir
Return the mesh sub-directory name (usually "polyMesh")
Definition: polyMesh.H:312
static unsigned int defaultPrecision()
Return the default precision.
Definition: IOstream.H:461
List< FixedList< label, 8 > > cells() const
Return the cells for filling the block.
Definition: blockCreate.C:349
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:256
bool optionFound(const word &opt) const
Return true if the named option is found.
Definition: argListI.H:108
static void noParallel()
Remove the parallel options.
Definition: argList.C:167
A bounding box defined in terms of the points at its extremities.
Definition: boundBox.H:58
const Cmpt & z() const
Definition: VectorI.H:87
Xfer< T > xferCopy(const T &)
Construct by copying the contents of the arg.
bool optionReadIfPresent(const word &opt, T &) const
Read a value from the named option if present.
Definition: argListI.H:198
const Cmpt & y() const
Definition: VectorI.H:81
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Definition: IOdictionary.H:52
virtual const pointField & points() const
Return raw points.
Definition: polyMesh.C:1003
iterator find(const Key &)
Find and return an iterator set at the hashedEntry.
Definition: HashTable.C:142
Functions used by OpenFOAM that are specific to POSIX compliant operating systems and need to be repl...
dynamicFvMesh & mesh
A multi-block mesh generator.
Definition: blockMesh.H:61
void removeFiles(const fileName &instanceDir) const
Remove all files from mesh instance.
Definition: polyMesh.C:1204
const pointZoneMesh & pointZones() const
Return point zone mesh.
Definition: polyMesh.H:460
bool isDir(const fileName &, const bool followLink=true)
Does the name exist as a DIRECTORY in the file system?
Definition: POSIX.C:528
A class for handling words, derived from string.
Definition: word.H:59
fileName dictPath
const cellZoneMesh & cellZones() const
Return cell zone mesh.
Definition: polyMesh.H:472
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:120
const word dictName("particleTrackDict")
void addZones(const List< pointZone *> &pz, const List< faceZone *> &fz, const List< cellZone *> &cz)
Add mesh zones.
Definition: polyMesh.C:946
const vectorField & cellCentres() const
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
const Cmpt & x() const
Definition: VectorI.H:75
void setSize(const label)
Reset size of PtrList. If extending the PtrList, new entries are.
Definition: PtrList.C:131
const word & system() const
Return system name.
Definition: TimePaths.H:114
Creates a single block of cells from point coordinates, numbers of cells in each direction and an exp...
Definition: block.H:63
bool rmDir(const fileName &)
Remove a directory and its contents.
Definition: POSIX.C:1032
static const char nl
Definition: Ostream.H:265
A subset of mesh cells.
Definition: cellZone.H:61
word defaultFacesName
Definition: readKivaGrid.H:460
void setInstance(const fileName &)
Set the instance for mesh files.
Definition: polyMeshIO.C:32
const fileName & instance() const
Definition: IOobject.H:392
label patchi
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:63
Direct mesh changes based on v1.3 polyTopoChange syntax.
label start() const
Return start label of this patch in the polyMesh face list.
Definition: polyPatch.H:300
messageStream Info
label nPoints() const
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:110
static void addNote(const string &)
Add extra notes for the usage information.
Definition: argList.C:151
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
const word & zoneName() const
Return the (optional) zone name.
Foam::argList args(argc, argv)
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
bool exists(const fileName &, const bool checkGzip=true, const bool followLink=true)
Does the name exist (as DIRECTORY or FILE) in the file system?
Definition: POSIX.C:509
void shrink()
Shrink the allocated table to approx. twice number of elements.
Definition: HashTable.C:500
Namespace for OpenFOAM.