FreeFOAM The Cross-Platform CFD Toolkit
hexBlock.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd.
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 Description
25 
26 \*---------------------------------------------------------------------------*/
27 
28 #include "hexBlock.H"
29 
30 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
31 
32 namespace Foam
33 {
34 
35 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
36 
37 label hexBlock::vtxLabel(label a, label b, label c) const
38 {
39  return (a + b*(xDim_ + 1) + c*(xDim_ + 1)*(yDim_ + 1));
40 }
41 
42 
43 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
44 
45 // Construct from components
46 hexBlock::hexBlock(const label nx, const label ny, const label nz)
47 :
48  xDim_(nx),
49  yDim_(ny),
50  zDim_(nz),
51  blockHandedness_(noPoints),
52  points_((xDim_ + 1)*(yDim_ + 1)*(zDim_ + 1))
53 {}
54 
55 
56 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
57 
58 void hexBlock::readPoints(Istream& is)
59 {
60  forAll (points_, i)
61  {
62  is >> points_[i].x() >> points_[i].y() >> points_[i].z();
63  }
64 
65  // Calculate the handedness of the block
66  vector i = points_[xDim_] - points_[0];
67  vector j = points_[(xDim_ + 1)*yDim_] - points_[0];
68  vector k = points_[(xDim_ + 1)*(yDim_ + 1)*zDim_] - points_[0];
69 
70  if (((i ^ j) & k) > 0)
71  {
72  Info << "right-handed block" << endl;
73  blockHandedness_ = right;
74  }
75  else
76  {
77  Info << "left-handed block" << endl;
78  blockHandedness_ = left;
79  }
80 }
81 
82 
84 {
85  labelListList result(xDim_*yDim_*zDim_);
86 
87  label cellNo = 0;
88 
89  if (blockHandedness_ == right)
90  {
91  for (label k = 0; k <= zDim_ - 1; k++)
92  {
93  for (label j = 0; j <= yDim_ - 1; j++)
94  {
95  for (label i = 0; i <= xDim_ - 1; i++)
96  {
97  labelList& hexLabels = result[cellNo];
98  hexLabels.setSize(8);
99 
100  hexLabels[0] = vtxLabel(i, j, k);
101  hexLabels[1] = vtxLabel(i+1, j, k);
102  hexLabels[2] = vtxLabel(i+1, j+1, k);
103  hexLabels[3] = vtxLabel(i, j+1, k);
104  hexLabels[4] = vtxLabel(i, j, k+1);
105  hexLabels[5] = vtxLabel(i+1, j, k+1);
106  hexLabels[6] = vtxLabel(i+1, j+1, k+1);
107  hexLabels[7] = vtxLabel(i, j+1, k+1);
108 
109  cellNo++;
110  }
111  }
112  }
113  }
114  else if (blockHandedness_ == left)
115  {
116  for (label k = 0; k <= zDim_ - 1; k++)
117  {
118  for (label j = 0; j <= yDim_ - 1; j++)
119  {
120  for (label i = 0; i <= xDim_ - 1; i++)
121  {
122  labelList& hexLabels = result[cellNo];
123  hexLabels.setSize(8);
124 
125  hexLabels[0] = vtxLabel(i, j, k+1);
126  hexLabels[1] = vtxLabel(i+1, j, k+1);
127  hexLabels[2] = vtxLabel(i+1, j+1, k+1);
128  hexLabels[3] = vtxLabel(i, j+1, k+1);
129  hexLabels[4] = vtxLabel(i, j, k);
130  hexLabels[5] = vtxLabel(i+1, j, k);
131  hexLabels[6] = vtxLabel(i+1, j+1, k);
132  hexLabels[7] = vtxLabel(i, j+1, k);
133 
134  cellNo++;
135  }
136  }
137  }
138  }
139  else
140  {
141  FatalErrorIn("hexBlock::cellShapes()")
142  << "Unable to determine block handedness as points "
143  << "have not been read in yet"
144  << abort(FatalError);
145  }
146 
147  return result;
148 }
149 
150 
151 // Return block patch faces given direction and range limits
152 // From the cfx manual: direction
153 // 0 = solid (3-D patch),
154 // 1 = high i, 2 = high j, 3 = high k
155 // 4 = low i, 5 = low j, 6 = low k
156 faceList hexBlock::patchFaces(const label direc, const labelList& range) const
157 {
158  if (range.size() != 6)
159  {
161  (
162  "patchFaces(const label direc, const labelList& range) const"
163  ) << "Invalid size of the range array: " << range.size()
164  << ". Should be 6 (xMin, xMax, yMin, yMax, zMin, zMax"
165  << abort(FatalError);
166  }
167 
168  label xMinRange = range[0];
169  label xMaxRange = range[1];
170  label yMinRange = range[2];
171  label yMaxRange = range[3];
172  label zMinRange = range[4];
173  label zMaxRange = range[5];
174 
175  faceList result(0);
176 
177  switch (direc)
178  {
179  case 1:
180  {
181  // high i = xmax
182 
183  result.setSize
184  (
185  (yMaxRange - yMinRange + 1)*(zMaxRange - zMinRange + 1)
186  );
187 
188  label p = 0;
189  for (label k = zMinRange - 1; k <= zMaxRange - 1; k++)
190  {
191  for (label j = yMinRange - 1; j <= yMaxRange - 1; j++)
192  {
193  result[p].setSize(4);
194 
195  // set the points
196  result[p][0] = vtxLabel(xDim_, j, k);
197  result[p][1] = vtxLabel(xDim_, j+1, k);
198  result[p][2] = vtxLabel(xDim_, j+1, k+1);
199  result[p][3] = vtxLabel(xDim_, j, k+1);
200 
201  p++;
202  }
203  }
204 
205  result.setSize(p);
206  break;
207  }
208 
209  case 2:
210  {
211  // high j = ymax
212  result.setSize
213  (
214  (xMaxRange - xMinRange + 1)*(zMaxRange - zMinRange + 1)
215  );
216 
217  label p = 0;
218  for (label i = xMinRange - 1; i <= xMaxRange - 1; i++)
219  {
220  for (label k = zMinRange - 1; k <= zMaxRange - 1; k++)
221  {
222  result[p].setSize(4);
223 
224  // set the points
225  result[p][0] = vtxLabel(i, yDim_, k);
226  result[p][1] = vtxLabel(i, yDim_, k + 1);
227  result[p][2] = vtxLabel(i + 1, yDim_, k + 1);
228  result[p][3] = vtxLabel(i + 1, yDim_, k);
229 
230  p++;
231  }
232  }
233 
234  result.setSize(p);
235  break;
236  }
237 
238  case 3:
239  {
240  // high k = zmax
241  result.setSize
242  (
243  (xMaxRange - xMinRange + 1)*(yMaxRange - yMinRange + 1)
244  );
245 
246  label p = 0;
247  for (label i = xMinRange - 1; i <= xMaxRange - 1; i++)
248  {
249  for (label j = yMinRange - 1; j <= yMaxRange - 1; j++)
250  {
251  result[p].setSize(4);
252 
253  // set the points
254  result[p][0] = vtxLabel(i, j, zDim_);
255  result[p][1] = vtxLabel(i + 1, j, zDim_);
256  result[p][2] = vtxLabel(i + 1, j + 1, zDim_);
257  result[p][3] = vtxLabel(i, j + 1, zDim_);
258 
259  p++;
260  }
261  }
262 
263  result.setSize(p);
264  break;
265  }
266 
267  case 4:
268  {
269  // low i = xmin
270  result.setSize
271  (
272  (yMaxRange - yMinRange + 1)*(zMaxRange - zMinRange + 1)
273  );
274 
275  label p = 0;
276  for (label k = zMinRange - 1; k <= zMaxRange - 1; k++)
277  {
278  for (label j = yMinRange - 1; j <= yMaxRange - 1; j++)
279  {
280  result[p].setSize(4);
281 
282  // set the points
283  result[p][0] = vtxLabel(0, j, k);
284  result[p][1] = vtxLabel(0, j, k + 1);
285  result[p][2] = vtxLabel(0, j + 1, k + 1);
286  result[p][3] = vtxLabel(0, j + 1, k);
287 
288  p++;
289  }
290  }
291 
292  result.setSize(p);
293  break;
294  }
295 
296  case 5:
297  {
298  // low j = ymin
299  result.setSize
300  (
301  (xMaxRange - xMinRange + 1)*(zMaxRange - zMinRange + 1)
302  );
303 
304  label p = 0;
305  for (label i = xMinRange - 1; i <= xMaxRange - 1; i++)
306  {
307  for (label k = zMinRange - 1; k <= zMaxRange - 1; k++)
308  {
309  result[p].setSize(4);
310 
311  // set the points
312  result[p][0] = vtxLabel(i, 0, k);
313  result[p][1] = vtxLabel(i + 1, 0, k);
314  result[p][2] = vtxLabel(i + 1, 0, k + 1);
315  result[p][3] = vtxLabel(i, 0, k + 1);
316 
317  p++;
318  }
319  }
320 
321  result.setSize(p);
322  break;
323  }
324 
325  case 6:
326  {
327  // low k = zmin
328  result.setSize
329  (
330  (xMaxRange - xMinRange + 1)*(yMaxRange - yMinRange + 1)
331  );
332 
333  label p = 0;
334  for (label i = xMinRange - 1; i <= xMaxRange - 1; i++)
335  {
336  for (label j = yMinRange - 1; j <= yMaxRange - 1; j++)
337  {
338  result[p].setSize(4);
339 
340  // set the points
341  result[p][0] = vtxLabel(i, j, 0);
342  result[p][1] = vtxLabel(i, j + 1, 0);
343  result[p][2] = vtxLabel(i + 1, j + 1, 0);
344  result[p][3] = vtxLabel(i + 1, j, 0);
345 
346  p++;
347  }
348  }
349 
350  result.setSize(p);
351  break;
352  }
353 
354  default:
355  {
357  (
358  "patchFaces(const label direc, const labelList& range) const"
359  ) << "direction out of range (1 to 6): " << direc
360  << abort(FatalError);
361  }
362  }
363 
364  // Correct the face orientation based on the handedness of the block.
365  // Do nothing for the right-handed block
366  if (blockHandedness_ == noPoints)
367  {
369  (
370  "patchFaces(const label direc, const labelList& range) const"
371  ) << "Unable to determine block handedness as points "
372  << "have not been read in yet"
373  << abort(FatalError);
374  }
375  else if (blockHandedness_ == left)
376  {
377  // turn all faces inside out
378  forAll (result, faceI)
379  {
380  result[faceI] = result[faceI].reverseFace();
381  }
382  }
383 
384  return result;
385 }
386 
387 
388 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
389 
390 } // End namespace Foam
391 
392 // ************************ vim: set sw=4 sts=4 et: ************************ //