PDRblock.H
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------* \
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2019-2020 OpenCFD Ltd.
9 -------------------------------------------------------------------------------
10 License
11  This file is part of OpenFOAM.
12 
13  OpenFOAM is free software: you can redistribute it and/or modify it
14  under the terms of the GNU General Public License as published by
15  the Free Software Foundation, either version 3 of the License, or
16  (at your option) any later version.
17 
18  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21  for more details.
22 
23  You should have received a copy of the GNU General Public License
24  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25 
26 Class
27  Foam::PDRblock
28 
29 Description
30  A single block x-y-z rectilinear mesh addressable as i,j,k with
31  simplified creation. Some of the input is similar to blockMeshDict,
32  but since this specialization is for a single-block that is aligned
33  with the x-y-z directions, it provides a different means of specifying
34  the mesh.
35 
36  Dictionary controls
37  \table
38  Property | Description | Required | Default
39  x | X-direction grid specification | yes |
40  y | Y-direction grid specification | yes |
41  z | Z-direction grid specification | yes |
42  scale | Point scaling | no | 1.0
43  expansion | Type of expansion (ratio/relative) | no | ratio
44  boundary | Boundary patches | yes |
45  defaultPatch | Default patch specification | no |
46  \endtable
47 
48  Grid coordinate controls
49  \table
50  Property| Description | Required | Default
51  points | Locations defining the mesh segment | yes |
52  nCells | Divisions per mesh segment | yes |
53  ratios | Expansion values per segment | no |
54  \endtable
55 
56  A negative expansion value is trapped and treated as its reciprocal.
57  by default, the expansion is as per blockMesh and represents the ratio
58  of end-size / start-size for the section.
59  Alternatively, the relative size can be given.
60 
61 SourceFiles
62  PDRblockI.H
63  PDRblock.C
64  PDRblockCreate.C
65 
66 \*---------------------------------------------------------------------------*/
67 
68 #ifndef PDRblock_H
69 #define PDRblock_H
70 
71 #include "ijkMesh.H"
72 #include "boundBox.H"
73 #include "pointField.H"
74 #include "faceList.H"
75 #include "Enum.H"
76 
77 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
78 
79 namespace Foam
80 {
81 
82 // Forward Declarations
83 class IOobject;
84 class polyMesh;
85 
86 /*---------------------------------------------------------------------------*\
87  Class PDRblock Declaration
88 \*---------------------------------------------------------------------------*/
89 
90 class PDRblock
91 :
92  public ijkMesh
93 {
94 public:
95 
96  // Data Types
97 
98  //- The expansion type
99  enum expansionType
100  {
102  EXPAND_RATIO,
104  };
105 
106 
107  // Public Classes
108 
109  //- Grid locations in an axis direction.
110  // The number of points is one larger than the number of elements
111  // it represents
112  class location
113  :
114  public scalarList
115  {
116  public:
117 
118  //- The locations are valid if they contain 2 or more points
119  inline bool valid() const;
120 
121  //- The number of cells in this direction.
122  inline label nCells() const;
123 
124  //- The number of points in this direction.
125  inline label nPoints() const;
126 
127  //- True if the location is within the range
128  inline bool contains(const scalar p) const;
129 
130  //- The first() value is considered the min value.
131  inline const scalar& min() const;
132 
133  //- The last() value is considered the max value.
134  inline const scalar& max() const;
135 
136  //- Mid-point location, zero for an empty list.
137  inline scalar centre() const;
138 
139  //- Check that element index is within valid range.
140  inline void checkIndex(const label i) const;
141 
142  //- Cell size at element position.
143  inline scalar width(const label i) const;
144 
145  //- Cell centre at element position.
146  // Treats -1 and nCells positions like a halo cell.
147  inline scalar C(const label i) const;
148 
149  //- Find the cell index enclosing this location
150  // \return -1 for out-of-bounds
151  label findCell(const scalar p) const;
152 
153  //- Find the grid index, within the given tolerance
154  // Return -1 for out-of-bounds and -2 for a point that is
155  // within bounds, but not aligned with a grid point.
156  label findIndex(const scalar p, const scalar tol) const;
157 
158  //- If out of range, return the respective min/max limits,
159  //- otherwise return the value itself.
160  // If the range is invalid, always return the value.
161  inline const scalar& clip(const scalar& val) const;
162  };
163 
164 
165 private:
166 
167  // Data Types
168 
169  //- Named enumerations for the expansion type
170  const static Enum<expansionType> expansionNames_;
171 
172 
173  // Private Classes
174 
175  //- Extracted patch settings
176  struct boundaryEntry
177  {
178  //- The patch name
179  word name_;
180 
181  //- The patch type
182  word type_;
183 
184  //- The patch size
185  label size_;
186 
187  //- The associated block face ids [0,5]
188  labelList faces_;
189  };
190 
191 
192  // Private Data
193 
194  //- Verbosity
195  bool verbose_;
196 
197  //- The grid points in all (i,j,k / x,y,z) directions.
198  Vector<location> grid_;
199 
200  //- The mesh bounding box
201  boundBox bounds_;
202 
203  //- The boundary patch information
204  PtrList<boundaryEntry> patches_;
205 
206  //- The min edge length
207  scalar minEdgeLen_;
208 
209 
210  // Private Member Functions
211 
212  //- Check that points increase monotonically
213  static bool checkMonotonic
214  (
215  const direction cmpt,
216  const UList<scalar>& pts
217  );
218 
219  //- Adjust sizing for updated grid points
220  void adjustSizes();
221 
222  //- Read and define grid points in given direction
223  void readGridControl
224  (
225  const direction cmpt,
226  const dictionary& dict,
227  const scalar scaleFactor = -1,
228  expansionType expandType = expansionType::EXPAND_RATIO
229  );
230 
231  //- Read "boundary" information
232  void readBoundary(const dictionary& dict);
233 
234  //- Populate point field for the block
235  void createPoints(pointField& pts) const;
236 
237  //- Add internal faces to lists.
238  // Lists must be properly sized!
239  // \return the number of faces added
240  label addInternalFaces
241  (
242  faceList::iterator& faceIter,
243  labelList::iterator& ownIter,
244  labelList::iterator& neiIter
245  ) const;
246 
247  //- Add boundary faces for the shape face to lists
248  // Lists must be properly sized!
249  // \return the number of faces added
250  label addBoundaryFaces
251  (
252  const direction shapeFacei,
253  faceList::iterator& faceIter,
254  labelList::iterator& ownIter
255  ) const;
256 
257  //- Obtain i,j,k index for cell enclosing this location
258  // \return false for out-of-bounds
259  bool findCell(const point& pt, labelVector& pos) const;
260 
261  //- Obtain i,j,k grid index for point location
262  // \return false for out-of-bounds and off-grid
263  bool gridIndex
264  (
265  const point& pt,
266  labelVector& pos,
267  const scalar tol
268  ) const;
269 
270 
271 public:
272 
273  // Static Member Functions
274 
275  //- Return a PDRblock reference to a nullObject
276  static const PDRblock& null();
277 
278 
279  // Constructors
280 
281  //- Construct zero-size
282  inline PDRblock();
283 
284  //- Construct from components
285  PDRblock
286  (
287  const UList<scalar>& xgrid,
288  const UList<scalar>& ygrid,
289  const UList<scalar>& zgrid
290  );
291 
292  //- Construct from dictionary
293  explicit PDRblock(const dictionary& dict, bool verbose=false);
294 
295 
296  // Member Functions
297 
298  //- Read dictionary
299  bool read(const dictionary& dict);
300 
301  //- Reset grid locations and mesh i-j-k sizing
302  void reset
303  (
304  const UList<scalar>& xgrid,
305  const UList<scalar>& ygrid,
306  const UList<scalar>& zgrid
307  );
308 
309 
310  // Access
311 
312  //- The grid point locations in the i,j,k (x,y,z) directions.
313  inline const Vector<location>& grid() const;
314 
315 
316  // Mesh information
317 
318  //- The mesh bounding box
319  inline const boundBox& bounds() const;
320 
321  //- The min edge length
322  inline const scalar& minEdgeLen() const;
323 
324  //- Cell size in x-direction at i position.
325  inline scalar dx(const label i) const;
326 
327  //- Cell size in x-direction at i position.
328  inline scalar dx(const labelVector& ijk) const;
329 
330  //- Cell size in y-direction at j position.
331  inline scalar dy(const label j) const;
332 
333  //- Cell size in y-direction at j position.
334  inline scalar dy(const labelVector& ijk) const;
335 
336  //- Cell size in z-direction at k position.
337  inline scalar dz(const label k) const;
338 
339  //- Cell size in z-direction at k position.
340  inline scalar dz(const labelVector& ijk) const;
341 
342  //- Cell dimensions at i,j,k position.
343  inline vector span(const label i, const label j, const label k) const;
344 
345  //- Cell dimensions at i,j,k position.
346  inline vector span(const labelVector& ijk) const;
347 
348  //- Grid point at i,j,k position.
349  inline point grid(const label i, const label j, const label k) const;
350 
351  //- Grid point at i,j,k position.
352  inline point grid(const labelVector& ijk) const;
353 
354  //- Cell centre at i,j,k position.
355  inline point C(const label i, const label j, const label k) const;
356 
357  //- Cell centre at i,j,k position.
358  inline point C(const labelVector& ijk) const;
359 
360  //- Cell volume at i,j,k position.
361  inline scalar V(const label i, const label j, const label k) const;
362 
363  //- Cell volume at i,j,k position.
364  inline scalar V(const labelVector& ijk) const;
365 
366  //- Characteristic cell size at i,j,k position.
367  // This is the cubic root of the volume
368  inline scalar width(const label i, const label j, const label k) const;
369 
370  //- Characteristic cell size at i,j,k position.
371  // This is the cubic root of the volume
372  inline scalar width(const labelVector& ijk) const;
373 
374 
375  // Searching
376 
377  //- Return i,j,k index for cell enclosing this location
378  // The value (-1,-1,-1) is returned for out-of-bounds (not found).
379  labelVector findCell(const point& pt) const;
380 
381  //- Obtain i,j,k grid index for point location within specified
382  // relative tolerance of the minEdgeLen.
383  // The value (-1,-1,-1) is returned for out-of-bounds (not found).
384  // and off-grid
385  labelVector gridIndex(const point& pt, const scalar relTol=0.01) const;
386 
387 
388  // Mesh generation
389 
390  //- Create polyMesh for grid definition and patch information
391  autoPtr<polyMesh> mesh(const IOobject& io) const;
392 
393 };
394 
395 
396 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
397 
398 } // End namespace Foam
399 
400 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
401 
402 #include "PDRblockI.H"
403 
404 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
405 
406 #endif
407 
408 // ************************************************************************* //
Foam::labelList
List< label > labelList
A List of labels.
Definition: List.H:71
Foam::PDRblock::width
scalar width(const label i, const label j, const label k) const
Characteristic cell size at i,j,k position.
Definition: PDRblockI.H:287
Foam::pointField
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:44
Foam::IOobject
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:104
Foam::PDRblock::dz
scalar dz(const label k) const
Cell size in z-direction at k position.
Definition: PDRblockI.H:194
Foam::scalarList
List< scalar > scalarList
A List of scalars.
Definition: scalarList.H:64
p
volScalarField & p
Definition: createFieldRefs.H:8
Foam::Enum< expansionType >
Foam::PDRblock::location::valid
bool valid() const
The locations are valid if they contain 2 or more points.
Definition: PDRblockI.H:43
Foam::PDRblock::EXPAND_RELATIVE
Relative expansion ratio.
Definition: PDRblock.H:162
Foam::word
A class for handling words, derived from Foam::string.
Definition: word.H:62
Foam::PDRblock::EXPAND_RATIO
End/start ratio.
Definition: PDRblock.H:161
Foam::PDRblock::dy
scalar dy(const label j) const
Cell size in y-direction at j position.
Definition: PDRblockI.H:182
PDRblockI.H
Foam::PDRblock::location::contains
bool contains(const scalar p) const
True if the location is within the range.
Definition: PDRblockI.H:61
Foam::PDRblock::minEdgeLen
const scalar & minEdgeLen() const
The min edge length.
Definition: PDRblockI.H:158
Foam::PDRblock::location::checkIndex
void checkIndex(const label i) const
Check that element index is within valid range.
Definition: PDRblockI.H:85
Foam::PDRblock::location::min
const scalar & min() const
The first() value is considered the min value.
Definition: PDRblockI.H:67
faceList.H
Foam::PDRblock::bounds
const boundBox & bounds() const
The mesh bounding box.
Definition: PDRblockI.H:164
Foam::PDRblock::expansionType
expansionType
The expansion type.
Definition: PDRblock.H:158
Foam::PDRblock::location::clip
const scalar & clip(const scalar &val) const
Definition: PDRblockI.H:133
Foam::labelVector
Vector< label > labelVector
Vector of labels.
Definition: labelVector.H:51
Foam::PDRblock::reset
void reset(const UList< scalar > &xgrid, const UList< scalar > &ygrid, const UList< scalar > &zgrid)
Reset grid locations and mesh i-j-k sizing.
Definition: PDRblock.C:583
Foam::PDRblock::location::findCell
label findCell(const scalar p) const
Find the cell index enclosing this location.
Definition: PDRblockLocation.C:33
Foam::PDRblock::mesh
autoPtr< polyMesh > mesh(const IOobject &io) const
Create polyMesh for grid definition and patch information.
Definition: PDRblockCreate.C:303
Foam::PDRblock::dx
scalar dx(const label i) const
Cell size in x-direction at i position.
Definition: PDRblockI.H:170
Foam::ijkMesh
A simple i-j-k (row-major order) to linear addressing for a rectilinear mesh. Since the underlying me...
Definition: ijkMesh.H:53
Foam::Field< vector >
ijkMesh.H
Foam::PDRblock::C
point C(const label i, const label j, const label k) const
Cell centre at i,j,k position.
Definition: PDRblockI.H:247
Foam::PtrList< boundaryEntry >
Foam::PDRblock::location::centre
scalar centre() const
Mid-point location, zero for an empty list.
Definition: PDRblockI.H:79
dict
dictionary dict
Definition: searchingEngine.H:14
Foam::dictionary
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:121
Foam::PDRblock::V
scalar V(const label i, const label j, const label k) const
Cell volume at i,j,k position.
Definition: PDRblockI.H:270
Foam::PDRblock::grid
const Vector< location > & grid() const
The grid point locations in the i,j,k (x,y,z) directions.
Definition: PDRblockI.H:152
Foam
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:33
Foam::PDRblock::location::nPoints
label nPoints() const
The number of points in this direction.
Definition: PDRblockI.H:55
Foam::vector
Vector< scalar > vector
A scalar version of the templated Vector.
Definition: vector.H:51
Foam::PDRblock
A single block x-y-z rectilinear mesh addressable as i,j,k with simplified creation....
Definition: PDRblock.H:149
Foam::PDRblock::location::width
scalar width(const label i) const
Cell size at element position.
Definition: PDRblockI.H:97
Foam::PDRblock::span
vector span(const label i, const label j, const label k) const
Cell dimensions at i,j,k position.
Definition: PDRblockI.H:207
pointField.H
boundBox.H
Foam::autoPtr
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
Definition: HashPtrTable.H:53
Foam::Vector
Templated 3D Vector derived from VectorSpace adding construction from 3 components,...
Definition: Vector.H:62
Foam::List< scalar >
Foam::PDRblock::location::findIndex
label findIndex(const scalar p, const scalar tol) const
Find the grid index, within the given tolerance.
Definition: PDRblockLocation.C:60
Foam::UList
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
Definition: HashTable.H:103
k
label k
Boltzmann constant.
Definition: LISASMDCalcMethod2.H:41
Foam::direction
uint8_t direction
Definition: direction.H:47
Foam::boundBox
A bounding box defined in terms of min/max extrema points.
Definition: boundBox.H:63
Foam::PDRblock::EXPAND_UNIFORM
Uniform expansion (ie, no expansion)
Definition: PDRblock.H:160
Foam::point
vector point
Point is a vector.
Definition: point.H:43
Foam::PDRblock::location::nCells
label nCells() const
The number of cells in this direction.
Definition: PDRblockI.H:49
Foam::PDRblock::read
bool read(const dictionary &dict)
Read dictionary.
Definition: PDRblock.C:556
Foam::PDRblock::location::max
const scalar & max() const
The last() value is considered the max value.
Definition: PDRblockI.H:73
Foam::PDRblock::PDRblock
PDRblock()
Construct zero-size.
Definition: PDRblockI.H:30
Foam::pos
dimensionedScalar pos(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:177
Enum.H
Foam::PDRblock::location::C
scalar C(const label i) const
Cell centre at element position.
Definition: PDRblockI.H:107