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-2025 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 "searchableSurfaceList.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  searchableSurfaceList 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,
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
154  (
155  const IOdictionary&,
156  const fileName& meshPath,
157  const word& regionName
158  );
159 
160  void check(const polyMesh&, const dictionary&) const;
161 
162  //- Determine the merge info and the final number of cells/points
163  void calcMergeInfo();
164 
165  //- Determine the merge info and the final number of cells/points
166  void calcMergeInfoFast();
167 
168  faceList createPatchFaces(const polyPatch& patchTopologyFaces) const;
169 
170  Pair<scalar> xCellSizes
171  (
172  const block& b,
173  const pointField& blockPoints,
174  const label j,
175  const label k
176  ) const;
177 
178  Pair<scalar> yCellSizes
179  (
180  const block& b,
181  const pointField& blockPoints,
182  const label i,
183  const label k
184  ) const;
185 
186  Pair<scalar> zCellSizes
187  (
188  const block& b,
189  const pointField& blockPoints,
190  const label i,
191  const label j
192  ) const;
193 
194  void printCellSizeRange(const Pair<scalar>& cellSizes) const;
195 
196  void printCellSizeRanges
197  (
198  const direction d,
199  const FixedList<Pair<scalar>, 4>& cellSizes
200  ) const;
201 
202  void createPoints() const;
203  void createCells() const;
204  void createPatches() const;
205 
206 
207 public:
208 
209  // Static Data Members
210 
211  ClassName("blockMesh");
212 
213  //- Switch checking block face orientation
214  // to ensure that all faces are outward-pointing.
215  // This check may fail if the block is intentionally very twisted
216  // for curved edges to be applied and can be switched off.
218 
219 
220  // Constructors
221 
222  //- Construct from IOdictionary, path and region name
223  blockMesh
224  (
225  const IOdictionary&,
226  const fileName& meshPath,
227  const word& regionName
228  );
229 
230  //- Disallow default bitwise copy construction
231  blockMesh(const blockMesh&) = delete;
232 
233 
234  //- Destructor
235  ~blockMesh();
236 
237 
238  // Member Functions
239 
240  // Access
241 
242  //- Access to input dictionary
243  const dictionary& meshDict() const
244  {
245  return meshDict_;
246  }
247 
248  //- Optional searchable geometry to project face-points to
249  const searchableSurfaceList& geometry() const
250  {
251  return geometry_;
252  }
253 
254  //- Reference to point field defining the blockMesh
255  // these points have not been scaled by scaleFactor
256  const pointField& vertices() const;
257 
258  //- Return the blockMesh topology as a polyMesh
259  const polyMesh& topology() const;
260 
261  //- Return the curved edges
262  const blockEdgeList& edges() const
263  {
264  return edges_;
265  }
266 
267  //- Return the curved faces
268  const blockFaceList& faces() const
269  {
270  return faces_;
271  }
272 
273  //- The scaling factor used to convert to meters
274  scalar scaleFactor() const;
275 
276  //- The points for the entire mesh
277  // these points have been scaled by scaleFactor
278  const pointField& points() const;
279 
280  //- Return cell shapes list
281  const cellShapeList& cells() const;
282 
283  //- Return the patch face lists
284  const faceListList& patches() const;
285 
286  //- Get patch information from the topology mesh
288 
289  //- Return patch names
290  wordList patchNames() const;
291 
292  //- Number of blocks with specified zones
293  label numZonedBlocks() const;
294 
295 
296  // Edit
297 
298  //- Enable/disable verbose information about the progress
299  void verbose(const bool on=true);
300 
301 
302  // Write
303 
304  //- Writes edges of blockMesh in OBJ format.
305  void writeTopology(Ostream&) const;
306 
307 
308  // Member Operators
309 
310  //- Disallow default bitwise assignment
311  void operator=(const blockMesh&) = delete;
312 };
313 
314 
315 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
316 
317 } // End namespace Foam
318 
319 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
320 
321 #endif
322 
323 // ************************************************************************* //
label k
A 1D vector of objects of type <T> with a fixed size <Size>.
Definition: FixedList.H:78
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Definition: IOdictionary.H:57
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:57
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
A simple wrapper around bool so that it can be read as a word: true/false, on/off,...
Definition: Switch.H:61
A multi-block mesh generator.
Definition: blockMesh.H:65
const polyMesh & topology() const
Return the blockMesh topology as a polyMesh.
const searchableSurfaceList & geometry() const
Optional searchable geometry to project face-points to.
Definition: blockMesh.H:248
ClassName("blockMesh")
void verbose(const bool on=true)
Enable/disable verbose information about the progress.
blockMesh(const IOdictionary &, const fileName &meshPath, const word &regionName)
Construct from IOdictionary, path and region name.
const dictionary & meshDict() const
Access to input dictionary.
Definition: blockMesh.H:242
static Switch checkBlockFaceOrientation
Switch checking block face orientation.
Definition: blockMesh.H:216
void operator=(const blockMesh &)=delete
Disallow default bitwise assignment.
label numZonedBlocks() const
Number of blocks with specified zones.
wordList patchNames() const
Return patch names.
const pointField & vertices() const
Reference to point field defining the blockMesh.
scalar scaleFactor() const
The scaling factor used to convert to meters.
PtrList< dictionary > patchDicts() const
Get patch information from the topology mesh.
const blockEdgeList & edges() const
Return the curved edges.
Definition: blockMesh.H:261
const faceListList & patches() const
Return the patch face lists.
~blockMesh()
Destructor.
const blockFaceList & faces() const
Return the curved faces.
Definition: blockMesh.H:267
void writeTopology(Ostream &) const
Writes edges of blockMesh in OBJ format.
const cellShapeList & cells() const
Return cell shapes list.
const pointField & points() const
The points for the entire mesh.
Creates a single block of cells from point coordinates, numbers of cells in each direction and an exp...
Definition: block.H:66
A list of keywords followed by any number of values (e.g. words and numbers) or sub-dictionaries.
Definition: dictionary.H:162
A class for handling file names.
Definition: fileName.H:82
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:80
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:70
Container for searchableSurfaces.
A class for handling words, derived from string.
Definition: word.H:62
volScalarField & b
Definition: createFields.H:27
Namespace for OpenFOAM.
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 & regionName(const solver &region)
Definition: solver.H:218
uint8_t direction
Definition: direction.H:45
wordList patchTypes(nPatches)
word meshPath
Definition: setMeshPath.H:1