blockMesh.H
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-2022 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 Class
25  Foam::blockMesh
26 
27 Description
28  A multi-block mesh generator
29 
30  Note:
31  The vertices, cells and patches for filling the blocks are
32  demand-driven.
33 
34 SourceFiles
35  blockMesh.C
36  blockMeshCheck.C
37  blockMeshCreate.C
38  blockMeshMerge.C
39  blockMeshTopology.C
40 
41 \*---------------------------------------------------------------------------*/
42 
43 #ifndef blockMesh_H
44 #define blockMesh_H
45 
46 #include "blockList.H"
47 #include "searchableSurfaces.H"
48 #include "polyMesh.H"
49 #include "IOdictionary.H"
50 #include "blockVertexList.H"
51 #include "blockEdgeList.H"
52 #include "blockFaceList.H"
53 
54 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
55 
56 namespace Foam
57 {
58 
59 /*---------------------------------------------------------------------------*\
60  Class blockMesh Declaration
61 \*---------------------------------------------------------------------------*/
62 
63 class blockMesh
64 :
65  public blockList
66 {
67  // Private Data
68 
69  //- Reference to mesh dictionary
70  const IOdictionary& meshDict_;
71 
72  //- Switch for verbose output
73  bool verboseOutput;
74 
75  //- Switch checking face consistency (defaults to true)
76  bool checkFaceCorrespondence_;
77 
78  //- Optional searchable geometry to project face-points to
79  searchableSurfaces geometry_;
80 
81  //- The scaling factor to convert to meters
82  scalar scaleFactor_;
83 
84  //- The list of block vertices
85  blockVertexList blockVertices_;
86 
87  //- The list of block vertex positions
88  pointField vertices_;
89 
90  //- The list of curved edges
91  blockEdgeList edges_;
92 
93  //- The list of curved faces
94  blockFaceList faces_;
95 
96  //- The blocks themselves (the topology) as a polyMesh
97  polyMesh* topologyPtr_;
98 
99  label nPoints_;
100 
101  //- The sum of all cells in each block
102  label nCells_;
103 
104  //- The point offset added to each block
105  labelList blockOffsets_;
106 
107  //- The merge points information
108  labelList mergeList_;
109 
110  mutable pointField points_;
111 
112  mutable cellShapeList cells_;
113 
114  mutable faceListList patches_;
115 
116 
117  // Private Member Functions
118 
119  template<class Source>
120  void checkPatchLabels
121  (
122  const Source& source,
123  const word& patchName,
124  const pointField& points,
125  faceList& patchShapes
126  ) const;
127 
128  void readPatches
129  (
130  const dictionary& meshDescription,
131  faceListList& tmpBlocksPatches,
134  wordList& nbrPatchNames
135  );
136 
137  void readBoundary
138  (
139  const dictionary& meshDescription,
140  wordList& patchNames,
141  faceListList& tmpBlocksPatches,
143  );
144 
145  void createCellShapes(cellShapeList& tmpBlockCells);
146 
147  void defaultPatchError
148  (
149  const word& defaultPatchName,
150  const dictionary& meshDescription
151  ) const;
152 
153  polyMesh* createTopology(const IOdictionary&, const word& regionName);
154 
155  void check(const polyMesh&, const dictionary&) const;
156 
157  //- Determine the merge info and the final number of cells/points
158  void calcMergeInfo();
159 
160  //- Determine the merge info and the final number of cells/points
161  void calcMergeInfoFast();
162 
163  faceList createPatchFaces(const polyPatch& patchTopologyFaces) const;
164 
165  Pair<scalar> xCellSizes
166  (
167  const block& b,
168  const pointField& blockPoints,
169  const label j,
170  const label k
171  ) const;
172 
173  Pair<scalar> yCellSizes
174  (
175  const block& b,
176  const pointField& blockPoints,
177  const label i,
178  const label k
179  ) const;
180 
181  Pair<scalar> zCellSizes
182  (
183  const block& b,
184  const pointField& blockPoints,
185  const label i,
186  const label j
187  ) const;
188 
189  void printCellSizeRange(const Pair<scalar>& cellSizes) const;
190 
191  void printCellSizeRanges
192  (
193  const direction d,
194  const FixedList<Pair<scalar>, 4>& cellSizes
195  ) const;
196 
197  void createPoints() const;
198  void createCells() const;
199  void createPatches() const;
200 
201 
202 public:
203 
204  // Static Data Members
205 
206  ClassName("blockMesh");
207 
208  //- Switch checking block face orientation
209  // to ensure that all faces are outward-pointing.
210  // This check may fail if the block is intentionally very twisted
211  // for curved edges to be applied and can be switched off.
213 
214 
215  // Constructors
216 
217  //- Construct from IOdictionary
218  blockMesh(const IOdictionary&, const word& regionName);
219 
220  //- Disallow default bitwise copy construction
221  blockMesh(const blockMesh&) = delete;
222 
223 
224  //- Destructor
225  ~blockMesh();
226 
227 
228  // Member Functions
229 
230  // Access
231 
232  //- Access to input dictionary
233  const dictionary& meshDict() const
234  {
235  return meshDict_;
236  }
237 
238  //- Optional searchable geometry to project face-points to
239  const searchableSurfaces& geometry() const
240  {
241  return geometry_;
242  }
243 
244  //- Reference to point field defining the blockMesh
245  // these points have not been scaled by scaleFactor
246  const pointField& vertices() const;
247 
248  //- Return the blockMesh topology as a polyMesh
249  const polyMesh& topology() const;
250 
251  //- Return the curved edges
252  const blockEdgeList& edges() const
253  {
254  return edges_;
255  }
256 
257  //- Return the curved faces
258  const blockFaceList& faces() const
259  {
260  return faces_;
261  }
262 
263  //- The scaling factor used to convert to meters
264  scalar scaleFactor() const;
265 
266  //- The points for the entire mesh
267  // these points have been scaled by scaleFactor
268  const pointField& points() const;
269 
270  //- Return cell shapes list
271  const cellShapeList& cells() const;
272 
273  //- Return the patch face lists
274  const faceListList& patches() const;
275 
276  //- Get patch information from the topology mesh
278 
279  //- Return patch names
280  wordList patchNames() const;
281 
282  //- Number of blocks with specified zones
283  label numZonedBlocks() const;
284 
285 
286  // Edit
287 
288  //- Enable/disable verbose information about the progress
289  void verbose(const bool on=true);
290 
291 
292  // Write
293 
294  //- Writes edges of blockMesh in OBJ format.
295  void writeTopology(Ostream&) const;
296 
297 
298  // Member Operators
299 
300  //- Disallow default bitwise assignment
301  void operator=(const blockMesh&) = delete;
302 };
303 
304 
305 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
306 
307 } // End namespace Foam
308 
309 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
310 
311 #endif
312 
313 // ************************************************************************* //
const blockFaceList & faces() const
Return the curved faces.
Definition: blockMesh.H:257
Foam::word regionName
A 1D vector of objects of type <T> with a fixed size <Size>.
Definition: FixedList.H:54
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
uint8_t direction
Definition: direction.H:45
const faceListList & patches() const
Return the patch face lists.
wordList patchTypes(nPatches)
const dimensionedScalar b
Wien displacement law constant: default SI units: [m K].
Definition: createFields.H:27
void writeTopology(Ostream &) const
Writes edges of blockMesh in OBJ format.
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
scalar scaleFactor() const
The scaling factor used to convert to meters.
const pointField & vertices() const
Reference to point field defining the blockMesh.
label k
Boltzmann constant.
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Definition: IOdictionary.H:53
blockMesh(const IOdictionary &, const word &regionName)
Construct from IOdictionary.
A multi-block mesh generator.
Definition: blockMesh.H:62
A class for handling words, derived from string.
Definition: word.H:59
const searchableSurfaces & geometry() const
Optional searchable geometry to project face-points to.
Definition: blockMesh.H:238
const polyMesh & topology() const
Return the blockMesh topology as a polyMesh.
Container for searchableSurfaces.
const blockEdgeList & edges() const
Return the curved edges.
Definition: blockMesh.H:251
void verbose(const bool on=true)
Enable/disable verbose information about the progress.
Creates a single block of cells from point coordinates, numbers of cells in each direction and an exp...
Definition: block.H:63
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:54
PtrList< dictionary > patchDicts() const
Get patch information from the topology mesh.
~blockMesh()
Destructor.
static Switch checkBlockFaceOrientation
Switch checking block face orientation.
Definition: blockMesh.H:211
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 cellShapeList & cells() const
Return cell shapes list.
const pointField & points() const
The points for the entire mesh.
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:76
label numZonedBlocks() const
Number of blocks with specified zones.
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
const dictionary & meshDict() const
Access to input dictionary.
Definition: blockMesh.H:232
wordList patchNames() const
Return patch names.
ClassName("blockMesh")
Namespace for OpenFOAM.
void operator=(const blockMesh &)=delete
Disallow default bitwise assignment.