blockMeshCreate.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-2019 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 \*---------------------------------------------------------------------------*/
25 
26 #include "blockMesh.H"
27 #include "cellModeller.H"
28 
29 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
30 
31 Foam::Pair<Foam::scalar> Foam::blockMesh::xCellSizes
32 (
33  const block& b,
34  const pointField& blockPoints,
35  const label j,
36  const label k
37 ) const
38 {
39  return Pair<scalar>
40  (
41  mag
42  (
43  blockPoints[b.pointLabel(1, j, k)]
44  - blockPoints[b.pointLabel(0, j, k)]
45  ),
46 
47  mag
48  (
49  blockPoints[b.pointLabel(b.density().x() - 1, j, k)]
50  - blockPoints[b.pointLabel(b.density().x(), j, k)]
51  )
52  );
53 }
54 
55 
56 Foam::Pair<Foam::scalar> Foam::blockMesh::yCellSizes
57 (
58  const block& b,
59  const pointField& blockPoints,
60  const label i,
61  const label k
62 ) const
63 {
64  return Pair<scalar>
65  (
66  mag
67  (
68  blockPoints[b.pointLabel(i, 0, k)]
69  - blockPoints[b.pointLabel(i, 1, k)]
70  ),
71 
72  mag
73  (
74  blockPoints[b.pointLabel(i, b.density().y() - 1, k)]
75  - blockPoints[b.pointLabel(i, b.density().y(), k)]
76  )
77  );
78 }
79 
80 
81 Foam::Pair<Foam::scalar> Foam::blockMesh::zCellSizes
82 (
83  const block& b,
84  const pointField& blockPoints,
85  const label i,
86  const label j
87 ) const
88 {
89  return Pair<scalar>
90  (
91  mag
92  (
93  blockPoints[b.pointLabel(i, j, 0)]
94  - blockPoints[b.pointLabel(i, j, 1)]
95  ),
96 
97  mag
98  (
99  blockPoints[b.pointLabel(i, j, b.density().z() - 1)]
100  - blockPoints[b.pointLabel(i, j, b.density().z())]
101  )
102  );
103 }
104 
105 
106 void Foam::blockMesh::printCellSizeRange
107 (
108  const Pair<scalar>& cellSizes
109 ) const
110 {
111  if (cellSizes.first() != cellSizes.second())
112  {
113  Info<< scaleFactor_*cellSizes.first() << " .. "
114  << scaleFactor_*cellSizes.second();
115  }
116  else
117  {
118  Info<< scaleFactor_*cellSizes.first();
119  }
120 }
121 
122 
123 void Foam::blockMesh::printCellSizeRanges
124 (
125  const direction d,
126  const FixedList<Pair<scalar>, 4>& cellSizes
127 ) const
128 {
129  static const char dNames[3] = {'i', 'j', 'k'};
130 
131  const scalar d0 = max
132  (
133  max(cellSizes[0].first(), cellSizes[0].second()),
134  small
135  );
136 
137  bool uniform = true;
138  forAll(cellSizes, i)
139  {
140  uniform = uniform
141  && mag(cellSizes[i].first() - cellSizes[0].first())/d0 < 1e-5
142  && mag(cellSizes[i].second() - cellSizes[0].second())/d0 < 1e-5;
143  }
144 
145  Info<< " " << dNames[d] << " : ";
146  if (uniform)
147  {
148  printCellSizeRange(cellSizes[0]);
149  }
150  else
151  {
152  forAll(cellSizes, i)
153  {
154  printCellSizeRange(cellSizes[i]);
155  Info<< " ";
156  }
157  }
158  Info<< endl;
159 }
160 
161 
162 void Foam::blockMesh::createPoints() const
163 {
164  const blockList& blocks = *this;
165 
166  if (verboseOutput)
167  {
168  Info<< "Creating points with scale " << scaleFactor_ << endl;
169  }
170 
171  points_.setSize(nPoints_);
172 
173  forAll(blocks, blocki)
174  {
175  const pointField& blockPoints = blocks[blocki].points();
176 
177  if (verboseOutput)
178  {
179  Info<< " Block " << blocki << " cell size :" << nl;
180 
181  printCellSizeRanges
182  (
183  0,
184  {
185  xCellSizes(blocks[blocki], blockPoints, 0, 0),
186  xCellSizes(blocks[blocki], blockPoints, 1, 0),
187  xCellSizes(blocks[blocki], blockPoints, 0, 1),
188  xCellSizes(blocks[blocki], blockPoints, 1, 1)
189  }
190  );
191 
192  printCellSizeRanges
193  (
194  1,
195  {
196  yCellSizes(blocks[blocki], blockPoints, 0, 0),
197  yCellSizes(blocks[blocki], blockPoints, 1, 0),
198  yCellSizes(blocks[blocki], blockPoints, 0, 1),
199  yCellSizes(blocks[blocki], blockPoints, 1, 1)
200  }
201  );
202 
203  printCellSizeRanges
204  (
205  2,
206  {
207  zCellSizes(blocks[blocki], blockPoints, 0, 0),
208  zCellSizes(blocks[blocki], blockPoints, 1, 0),
209  zCellSizes(blocks[blocki], blockPoints, 0, 1),
210  zCellSizes(blocks[blocki], blockPoints, 1, 1)
211  }
212  );
213  }
214 
215  forAll(blockPoints, blockPointi)
216  {
217  points_
218  [
219  mergeList_
220  [
221  blockOffsets_[blocki] + blockPointi
222  ]
223  ] = scaleFactor_ * blockPoints[blockPointi];
224  }
225  }
226 }
227 
228 
229 void Foam::blockMesh::createCells() const
230 {
231  const blockList& blocks = *this;
232  const cellModel& hex = *(cellModeller::lookup("hex"));
233 
234  if (verboseOutput)
235  {
236  Info<< "Creating cells" << endl;
237  }
238 
239  cells_.setSize(nCells_);
240 
241  label cellLabel = 0;
242 
243  forAll(blocks, blocki)
244  {
245  const List<FixedList<label, 8>> blockCells(blocks[blocki].cells());
246 
247  forAll(blockCells, blockCelli)
248  {
249  labelList cellPoints(blockCells[blockCelli].size());
250 
251  forAll(cellPoints, cellPointi)
252  {
253  cellPoints[cellPointi] =
254  mergeList_
255  [
256  blockCells[blockCelli][cellPointi]
257  + blockOffsets_[blocki]
258  ];
259  }
260 
261  // Construct collapsed cell and add to list
262  cells_[cellLabel] = cellShape(hex, cellPoints, true);
263 
264  cellLabel++;
265  }
266  }
267 }
268 
269 
270 Foam::faceList Foam::blockMesh::createPatchFaces
271 (
272  const polyPatch& patchTopologyFaces
273 ) const
274 {
275  const blockList& blocks = *this;
276 
277  labelList blockLabels = patchTopologyFaces.polyPatch::faceCells();
278 
279  label nFaces = 0;
280 
281  forAll(patchTopologyFaces, patchTopologyFaceLabel)
282  {
283  const label blocki = blockLabels[patchTopologyFaceLabel];
284 
285  faceList blockFaces = blocks[blocki].blockShape().faces();
286 
287  forAll(blockFaces, blockFaceLabel)
288  {
289  if
290  (
291  blockFaces[blockFaceLabel]
292  == patchTopologyFaces[patchTopologyFaceLabel]
293  )
294  {
295  nFaces +=
296  blocks[blocki].boundaryPatches()[blockFaceLabel].size();
297  }
298  }
299  }
300 
301 
302  faceList patchFaces(nFaces);
303  face quadFace(4);
304  label faceLabel = 0;
305 
306  forAll(patchTopologyFaces, patchTopologyFaceLabel)
307  {
308  const label blocki = blockLabels[patchTopologyFaceLabel];
309 
310  faceList blockFaces = blocks[blocki].blockShape().faces();
311 
312  forAll(blockFaces, blockFaceLabel)
313  {
314  if
315  (
316  blockFaces[blockFaceLabel]
317  == patchTopologyFaces[patchTopologyFaceLabel]
318  )
319  {
320  const List<FixedList<label, 4>>& blockPatchFaces =
321  blocks[blocki].boundaryPatches()[blockFaceLabel];
322 
323  forAll(blockPatchFaces, blockFaceLabel)
324  {
325  // Lookup the face points
326  // and collapse duplicate point labels
327 
328  quadFace[0] =
329  mergeList_
330  [
331  blockPatchFaces[blockFaceLabel][0]
332  + blockOffsets_[blocki]
333  ];
334 
335  label nUnique = 1;
336 
337  for
338  (
339  label facePointLabel = 1;
340  facePointLabel < 4;
341  facePointLabel++
342  )
343  {
344  quadFace[nUnique] =
345  mergeList_
346  [
347  blockPatchFaces[blockFaceLabel][facePointLabel]
348  + blockOffsets_[blocki]
349  ];
350 
351  if (quadFace[nUnique] != quadFace[nUnique-1])
352  {
353  nUnique++;
354  }
355  }
356 
357  if (quadFace[nUnique-1] == quadFace[0])
358  {
359  nUnique--;
360  }
361 
362  if (nUnique == 4)
363  {
364  patchFaces[faceLabel++] = quadFace;
365  }
366  else if (nUnique == 3)
367  {
368  patchFaces[faceLabel++] = face
369  (
371  );
372  }
373  // else the face has collapsed to an edge or point
374  }
375  }
376  }
377  }
378 
379  patchFaces.setSize(faceLabel);
380 
381  return patchFaces;
382 }
383 
384 
385 void Foam::blockMesh::createPatches() const
386 {
387  const polyPatchList& topoPatches = topology().boundaryMesh();
388 
389  if (verboseOutput)
390  {
391  Info<< "Creating patches" << endl;
392  }
393 
394  patches_.setSize(topoPatches.size());
395 
396  forAll(topoPatches, patchi)
397  {
398  patches_[patchi] = createPatchFaces(topoPatches[patchi]);
399  }
400 }
401 
402 
403 // ************************************************************************* //
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:434
PtrList< polyPatch > polyPatchList
container classes for polyPatch
Definition: polyPatchList.H:45
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
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
static const cellModel * lookup(const word &)
Look up a model by name and return a pointer to the model or nullptr.
Definition: cellModeller.C:100
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: HashTable.H:59
uint8_t direction
Definition: direction.H:45
List< face > faceList
Definition: faceListFwd.H:43
PtrList< block > blockList
A PtrList of blocks.
Definition: blockList.H:45
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
SubList< label > subList
Declare type of subList.
Definition: List.H:199
label k
Boltzmann constant.
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:42
An ordered pair of two objects of type <T> with first() and second() elements.
Definition: contiguous.H:49
const polyMesh & topology() const
Return the blockMesh topology as a polyMesh.
List< label > labelList
A List of labels.
Definition: labelList.H:56
static const char nl
Definition: Ostream.H:260
face quadFace(4)
label size() const
Return the number of elements in the UPtrList.
Definition: UPtrListI.H:29
void setSize(const label)
Reset size of List.
Definition: List.C:281
label patchi
const cellShapeList & cells() const
Return cell shapes list.
messageStream Info
dimensioned< scalar > mag(const dimensioned< Type > &)
const doubleScalar e
Elementary charge.
Definition: doubleScalar.H:105
T & first()
Return reference to the first element of the list.
Definition: UPtrListI.H:43