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-2024 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 "argList.H"
52 #include "Time.H"
53 #include "systemDict.H"
54 #include "blockMesh.H"
55 #include "mergePatchPairs.H"
56 #include "polyTopoChange.H"
57 #include "emptyPolyPatch.H"
58 #include "cyclicPolyPatch.H"
59 #include "OFstream.H"
60 
61 using namespace Foam;
62 
63 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
64 
65 int main(int argc, char *argv[])
66 {
68  Foam::argList::removeOption("noFunctionObjects");
69  #include "addDictOption.H"
71  (
72  "blockTopology",
73  "write block edges and centres as .obj files"
74  );
76  (
77  "noClean",
78  "keep the existing files in the polyMesh"
79  );
80 
82  (
83  "Block description\n"
84  "\n"
85  " For a given block, the correspondence between the ordering of\n"
86  " vertex labels and face labels is shown below.\n"
87  " For vertex numbering in the sequence 0 to 7 (block, centre):\n"
88  " faces 0 (f0) and 1 are left and right, respectively;\n"
89  " faces 2 and 3 are front and back; \n"
90  " and faces 4 and 5 are bottom and top::\n"
91  "\n"
92  " 7 ---- 6\n"
93  " f5 |\\ |\\ f3\n"
94  " | | 4 ---- 5 \\\n"
95  " | 3 |--- 2 | \\\n"
96  " | \\| \\| f2\n"
97  " f4 0 ---- 1\n"
98  "\n"
99  " Z f0 ----- f1\n"
100  " | Y\n"
101  " | /\n"
102  " O --- X\n"
103  );
104 
105  #include "addRegionOption.H"
106  #include "setRootCase.H"
107  #include "createTime.H"
108 
109  const word dictName("blockMeshDict");
110 
112  word regionPath;
113 
114  // Check if the region is specified otherwise mesh the default region
116  {
117  Info<< nl << "Generating mesh for region " << regionName << endl;
118  regionPath = regionName;
119  }
120 
121  if (!args.optionFound("noClean"))
122  {
123  fileName polyMeshPath
124  (
125  runTime.path()/runTime.constant()/regionPath/polyMesh::meshSubDir
126  );
127 
128  if (exists(polyMeshPath))
129  {
130  if (exists(polyMeshPath/dictName))
131  {
132  Info<< "Not deleting polyMesh directory " << nl
133  << " " << polyMeshPath << nl
134  << " because it contains " << dictName << endl;
135  }
136  else
137  {
138  Info<< "Deleting polyMesh directory" << nl
139  << " " << polyMeshPath << endl;
140  rmDir(polyMeshPath);
141  }
142  }
143  }
144 
145  typeIOobject<IOdictionary> meshDictIO
146  (
147  systemDictIO(dictName, args, runTime, regionName)
148  );
149 
150  if (!meshDictIO.headerOk())
151  {
153  << "Cannot find file " << meshDictIO.relativeObjectPath()
154  << nl
155  << exit(FatalError);
156  }
157 
158  Info<< "Creating block mesh from\n "
159  << meshDictIO.relativeObjectPath() << endl;
160 
161  IOdictionary meshDict(meshDictIO);
162  blockMesh blocks(meshDict, regionName);
163 
164 
165  if (args.optionFound("blockTopology"))
166  {
167  // Write mesh as edges.
168  {
169  fileName objMeshFile("blockTopology.obj");
170 
171  OFstream str(runTime.path()/objMeshFile);
172 
173  Info<< nl << "Dumping block structure as Lightwave obj format"
174  << " to " << objMeshFile << endl;
175 
176  blocks.writeTopology(str);
177  }
178 
179  // Write centres of blocks
180  {
181  fileName objCcFile("blockCentres.obj");
182 
183  OFstream str(runTime.path()/objCcFile);
184 
185  Info<< nl << "Dumping block centres as Lightwave obj format"
186  << " to " << objCcFile << endl;
187 
188  const polyMesh& topo = blocks.topology();
189 
190  const pointField& cellCentres = topo.cellCentres();
191 
192  forAll(cellCentres, celli)
193  {
194  // point cc = b.blockShape().centre(b.points());
195  const point& cc = cellCentres[celli];
196 
197  str << "v " << cc.x() << ' ' << cc.y() << ' ' << cc.z() << nl;
198  }
199  }
200 
201  Info<< nl << "end" << endl;
202 
203  return 0;
204  }
205 
206 
207  Info<< nl << "Creating polyMesh from blockMesh" << endl;
208 
209  word defaultFacesName = "defaultFaces";
210  word defaultFacesType = emptyPolyPatch::typeName;
211  polyMesh mesh
212  (
213  IOobject
214  (
215  regionName,
216  runTime.constant(),
217  runTime
218  ),
219  clone(blocks.points()), // could we reuse space?
220  blocks.cells(),
221  blocks.patches(),
222  blocks.patchNames(),
223  blocks.patchDicts(),
226  );
227 
228 
229  // Read in a list of dictionaries for the merge patch pairs
230  if (meshDict.found("mergePatchPairs"))
231  {
232  List<Pair<word>> patchPairNames
233  (
234  meshDict.lookup("mergePatchPairs")
235  );
236 
237  if (patchPairNames.size())
238  {
239  const word oldInstance = mesh.pointsInstance();
240  mergePatchPairs(mesh, patchPairNames, 1e-4);
241  mesh.setInstance(oldInstance);
242  }
243  }
244  else
245  {
246  Info<< nl << "No patch pairs to merge" << endl;
247  }
248 
249  label nZones = blocks.numZonedBlocks();
250 
251  if (nZones > 0)
252  {
253  Info<< nl << "Adding cell zones" << endl;
254 
255  // Map from zoneName to cellZone index
256  HashTable<label> zoneMap(nZones);
257 
258  // Cells per zone.
259  List<DynamicList<label>> zoneCells(nZones);
260 
261  // Running cell counter
262  label celli = 0;
263 
264  // Largest zone so far
265  label freeZoneI = 0;
266 
267  forAll(blocks, blockI)
268  {
269  const block& b = blocks[blockI];
270  const List<FixedList<label, 8>> blockCells = b.cells();
271  const word& zoneName = b.zoneName();
272 
273  if (zoneName.size())
274  {
275  HashTable<label>::const_iterator iter = zoneMap.find(zoneName);
276 
277  label zoneI;
278 
279  if (iter == zoneMap.end())
280  {
281  zoneI = freeZoneI++;
282 
283  Info<< " " << zoneI << '\t' << zoneName << endl;
284 
285  zoneMap.insert(zoneName, zoneI);
286  }
287  else
288  {
289  zoneI = iter();
290  }
291 
292  forAll(blockCells, i)
293  {
294  zoneCells[zoneI].append(celli++);
295  }
296  }
297  else
298  {
299  celli += blockCells.size();
300  }
301  }
302 
303 
304  List<cellZone*> cz(zoneMap.size());
305 
306  forAllConstIter(HashTable<label>, zoneMap, iter)
307  {
308  label zoneI = iter();
309 
310  cz[zoneI] = new cellZone
311  (
312  iter.key(),
313  zoneCells[zoneI].shrink(),
314  mesh.cellZones()
315  );
316  }
317 
318  mesh.pointZones().setSize(0);
319  mesh.faceZones().setSize(0);
320  mesh.cellZones().setSize(0);
321  mesh.addZones(List<pointZone*>(0), List<faceZone*>(0), cz);
322  }
323 
324 
325  // Detect any cyclic patches and force re-ordering of the faces
326  {
327  const polyPatchList& patches = mesh.boundaryMesh();
328  bool hasCyclic = false;
330  {
331  if (isA<cyclicPolyPatch>(patches[patchi]))
332  {
333  hasCyclic = true;
334  break;
335  }
336  }
337 
338  if (hasCyclic)
339  {
340  Info<< nl << "Detected cyclic patches; ordering boundary faces"
341  << endl;
342  const word oldInstance = mesh.instance();
343  polyTopoChange meshMod(mesh);
344  meshMod.changeMesh(mesh);
345  mesh.setInstance(oldInstance);
346  }
347  }
348 
349 
350  // Set the precision of the points data to 10
352 
353  Info<< nl << "Writing polyMesh" << endl;
354  mesh.removeFiles();
355  if (!mesh.write())
356  {
358  << "Failed writing polyMesh."
359  << exit(FatalError);
360  }
361 
362 
363  // Write summary
364  {
365  const polyPatchList& patches = mesh.boundaryMesh();
366 
367  Info<< "----------------" << nl
368  << "Mesh Information" << nl
369  << "----------------" << nl
370  << " " << "boundingBox: " << boundBox(mesh.points()) << nl
371  << " " << "nPoints: " << mesh.nPoints() << nl
372  << " " << "nCells: " << mesh.nCells() << nl
373  << " " << "nFaces: " << mesh.nFaces() << nl
374  << " " << "nInternalFaces: " << mesh.nInternalFaces() << nl;
375 
376  Info<< "----------------" << nl
377  << "Patches" << nl
378  << "----------------" << nl;
379 
381  {
382  const polyPatch& p = patches[patchi];
383 
384  Info<< " " << "patch " << patchi
385  << " (start: " << p.start()
386  << " size: " << p.size()
387  << ") name: " << p.name()
388  << nl;
389  }
390  }
391 
392  Info<< "\nEnd\n" << endl;
393 
394  return 0;
395 }
396 
397 
398 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
#define forAllConstIter(Container, container, iter)
Iterate across all elements in the container object of type.
Definition: UList.H:477
void shrink()
Shrink the allocated table to approx. twice number of elements.
Definition: HashTable.C:574
iterator find(const Key &)
Find and return an iterator set at the hashedEntry.
Definition: HashTable.C:167
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Definition: IOdictionary.H:57
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:99
const word & name() const
Return name.
Definition: IOobject.H:310
static unsigned int defaultPrecision()
Return the default precision.
Definition: IOstream.H:458
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: List.H:91
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
Output to file stream.
Definition: OFstream.H:86
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: PtrList.H:75
const Cmpt & z() const
Definition: VectorI.H:87
const Cmpt & y() const
Definition: VectorI.H:81
const Cmpt & x() const
Definition: VectorI.H:75
static void addNote(const string &)
Add extra notes for the usage information.
Definition: argList.C:159
static void addBoolOption(const word &opt, const string &usage="")
Add to a bool option to validOptions with usage information.
Definition: argList.C:118
bool optionFound(const word &opt) const
Return true if the named option is found.
Definition: argListI.H:114
bool optionReadIfPresent(const word &opt, T &) const
Read a value from the named option if present.
Definition: argListI.H:204
static void noParallel()
Remove the parallel options.
Definition: argList.C:175
static void removeOption(const word &opt)
Remove option from validOptions and from optionUsage.
Definition: argList.C:168
A multi-block mesh generator.
Definition: blockMesh.H:65
Creates a single block of cells from point coordinates, numbers of cells in each direction and an exp...
Definition: block.H:66
A bounding box defined in terms of the points at its extremities.
Definition: boundBox.H:59
A subset of mesh cells.
Definition: cellZone.H:58
A class for handling file names.
Definition: fileName.H:82
Class to stitch mesh by merging patch-pairs.
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:80
static word defaultRegion
Return the default region name.
Definition: polyMesh.H:267
static word meshSubDir
Return the mesh sub-directory name (usually "polyMesh")
Definition: polyMesh.H:270
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:70
Direct mesh changes based on v1.3 polyTopoChange syntax.
const vectorField & cellCentres() const
Templated form of IOobject providing type information for file reading and header type checking.
Definition: IOobject.H:531
A class for handling words, derived from string.
Definition: word.H:62
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:334
int main(int argc, char *argv[])
Definition: financialFoam.C:44
label patchi
const fvPatchList & patches
volScalarField & b
Definition: createFields.H:25
Namespace for OpenFOAM.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
const doubleScalar e
Definition: doubleScalar.H:106
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
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:257
const word & regionName(const solver &region)
Definition: solver.H:209
bool exists(const fileName &, const bool checkVariants=true, const bool followLink=true)
Does the name exist (as directory or file) in the file system?
Definition: POSIX.C:520
IOobject systemDictIO(const word &dictName, const argList &args, const objectRegistry &ob, const word &regionName=polyMesh::defaultRegion)
Definition: systemDict.C:33
messageStream Info
bool rmDir(const fileName &)
Remove a directory and its contents.
Definition: POSIX.C:1047
T clone(const T &t)
Definition: List.H:55
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
error FatalError
static const char nl
Definition: Ostream.H:266
word defaultFacesName
Definition: readKivaGrid.H:454
word defaultFacesType
Definition: readKivaGrid.H:455
Foam::argList args(argc, argv)
volScalarField & p