FreeFOAM The Cross-Platform CFD Toolkit
extrudedMesh.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 \*---------------------------------------------------------------------------*/
25 
26 #include "extrudedMesh.H"
27 #include <OpenFOAM/wallPolyPatch.H>
28 #include <meshTools/meshTools.H>
29 #include <OpenFOAM/ListOps.H>
30 #include <OpenFOAM/OFstream.H>
31 
32 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
33 
34 bool Foam::extrudedMesh::sameOrder(const face& f, const edge& e)
35 {
36  label i = findIndex(f, e[0]);
37 
38  label nextI = (i == f.size()-1 ? 0 : i+1);
39 
40  return f[nextI] == e[1];
41 }
42 
43 
44 template
45 <
46  class Face,
47  template<class> class FaceList,
48  class PointField
49 >
50 Foam::Xfer<Foam::pointField> Foam::extrudedMesh::extrudedPoints
51 (
52  const PrimitivePatch<Face, FaceList, PointField>& extrudePatch,
53  const extrudeModel& model
54 )
55 {
56  const pointField& surfacePoints = extrudePatch.localPoints();
57  const vectorField& surfaceNormals = extrudePatch.pointNormals();
58 
59  const label nLayers = model.nLayers();
60 
61  pointField ePoints((nLayers + 1)*surfacePoints.size());
62 
63  for (label layer=0; layer<=nLayers; layer++)
64  {
65  label offset = layer*surfacePoints.size();
66 
67  forAll(surfacePoints, i)
68  {
69  ePoints[offset + i] = model
70  (
71  surfacePoints[i],
72  surfaceNormals[i],
73  layer
74  );
75  }
76  }
77 
78  // return points for transferring
79  return xferMove(ePoints);
80 }
81 
82 
83 template<class Face, template<class> class FaceList, class PointField>
84 Foam::Xfer<Foam::faceList> Foam::extrudedMesh::extrudedFaces
85 (
86  const PrimitivePatch<Face, FaceList, PointField>& extrudePatch,
87  const extrudeModel& model
88 )
89 {
90  const pointField& surfacePoints = extrudePatch.localPoints();
91  const List<face>& surfaceFaces = extrudePatch.localFaces();
92  const edgeList& surfaceEdges = extrudePatch.edges();
93  const label nInternalEdges = extrudePatch.nInternalEdges();
94 
95  const label nLayers = model.nLayers();
96 
97  label nFaces =
98  (nLayers + 1)*surfaceFaces.size() + nLayers*surfaceEdges.size();
99 
100  faceList eFaces(nFaces);
101 
102  labelList quad(4);
103  label facei = 0;
104 
105  // Internal faces
106  for (label layer=0; layer<nLayers; layer++)
107  {
108  label currentLayerOffset = layer*surfacePoints.size();
109  label nextLayerOffset = currentLayerOffset + surfacePoints.size();
110 
111  // Vertical faces from layer to layer+1
112  for (label edgeI=0; edgeI<nInternalEdges; edgeI++)
113  {
114  const edge& e = surfaceEdges[edgeI];
115  const labelList& edgeFaces = extrudePatch.edgeFaces()[edgeI];
116 
117  face& f = eFaces[facei++];
118  f.setSize(4);
119 
120  if
121  (
122  (edgeFaces[0] < edgeFaces[1])
123  == sameOrder(surfaceFaces[edgeFaces[0]], e)
124  )
125  {
126  f[0] = e[0] + currentLayerOffset;
127  f[1] = e[1] + currentLayerOffset;
128  f[2] = e[1] + nextLayerOffset;
129  f[3] = e[0] + nextLayerOffset;
130  }
131  else
132  {
133  f[0] = e[1] + currentLayerOffset;
134  f[1] = e[0] + currentLayerOffset;
135  f[2] = e[0] + nextLayerOffset;
136  f[3] = e[1] + nextLayerOffset;
137  }
138  }
139 
140  // Faces between layer and layer+1
141  if (layer < nLayers-1)
142  {
143  forAll(surfaceFaces, i)
144  {
145  eFaces[facei++] =
146  face
147  (
148  surfaceFaces[i] //.reverseFace()
149  + nextLayerOffset
150  );
151  }
152  }
153  }
154 
155  // External side faces
156  for (label layer=0; layer<nLayers; layer++)
157  {
158  label currentLayerOffset = layer*surfacePoints.size();
159  label nextLayerOffset = currentLayerOffset + surfacePoints.size();
160 
161  // Side faces across layer
162  for (label edgeI=nInternalEdges; edgeI<surfaceEdges.size(); edgeI++)
163  {
164  const edge& e = surfaceEdges[edgeI];
165  const labelList& edgeFaces = extrudePatch.edgeFaces()[edgeI];
166 
167  face& f = eFaces[facei++];
168  f.setSize(4);
169 
170  if (sameOrder(surfaceFaces[edgeFaces[0]], e))
171  {
172  f[0] = e[0] + currentLayerOffset;
173  f[1] = e[1] + currentLayerOffset;
174  f[2] = e[1] + nextLayerOffset;
175  f[3] = e[0] + nextLayerOffset;
176  }
177  else
178  {
179  f[0] = e[1] + currentLayerOffset;
180  f[1] = e[0] + currentLayerOffset;
181  f[2] = e[0] + nextLayerOffset;
182  f[3] = e[1] + nextLayerOffset;
183  }
184  }
185  }
186 
187  // Bottom faces
188  forAll(surfaceFaces, i)
189  {
190  eFaces[facei++] = face(surfaceFaces[i]).reverseFace();
191  }
192 
193  // Top faces
194  forAll(surfaceFaces, i)
195  {
196  eFaces[facei++] =
197  face
198  (
199  surfaceFaces[i]
200  + nLayers*surfacePoints.size()
201  );
202  }
203 
204  // return points for transferring
205  return xferMove(eFaces);
206 }
207 
208 
209 template<class Face, template<class> class FaceList, class PointField>
210 Foam::Xfer<Foam::cellList> Foam::extrudedMesh::extrudedCells
211 (
212  const PrimitivePatch<Face, FaceList, PointField>& extrudePatch,
213  const extrudeModel& model
214 )
215 {
216  const List<face>& surfaceFaces = extrudePatch.localFaces();
217  const edgeList& surfaceEdges = extrudePatch.edges();
218  const label nInternalEdges = extrudePatch.nInternalEdges();
219 
220  const label nLayers = model.nLayers();
221 
222  cellList eCells(nLayers*surfaceFaces.size());
223 
224  // Size the cells
225  forAll(surfaceFaces, i)
226  {
227  const face& f = surfaceFaces[i];
228 
229  for (label layer=0; layer<nLayers; layer++)
230  {
231  eCells[i + layer*surfaceFaces.size()].setSize(f.size() + 2);
232  }
233  }
234 
235  // Current face count per cell.
236  labelList nCellFaces(eCells.size(), 0);
237 
238 
239  label facei = 0;
240 
241  for (label layer=0; layer<nLayers; layer++)
242  {
243  // Side faces from layer to layer+1
244  for (label i=0; i<nInternalEdges; i++)
245  {
246  // Get patch faces using edge
247  const labelList& edgeFaces = extrudePatch.edgeFaces()[i];
248 
249  // Get cells on this layer
250  label cell0 = layer*surfaceFaces.size() + edgeFaces[0];
251  label cell1 = layer*surfaceFaces.size() + edgeFaces[1];
252 
253  eCells[cell0][nCellFaces[cell0]++] = facei;
254  eCells[cell1][nCellFaces[cell1]++] = facei;
255 
256  facei++;
257  }
258 
259  // Faces between layer and layer+1
260  if (layer < nLayers-1)
261  {
262  forAll(surfaceFaces, i)
263  {
264  label cell0 = layer*surfaceFaces.size() + i;
265  label cell1 = (layer+1)*surfaceFaces.size() + i;
266 
267  eCells[cell0][nCellFaces[cell0]++] = facei;
268  eCells[cell1][nCellFaces[cell1]++] = facei;
269 
270  facei++;
271  }
272  }
273  }
274 
275  // External side faces
276  for (label layer=0; layer<nLayers; layer++)
277  {
278  // Side faces across layer
279  for (label i=nInternalEdges; i<surfaceEdges.size(); i++)
280  {
281  // Get patch faces using edge
282  const labelList& edgeFaces = extrudePatch.edgeFaces()[i];
283 
284  // Get cells on this layer
285  label cell0 = layer*surfaceFaces.size() + edgeFaces[0];
286 
287  eCells[cell0][nCellFaces[cell0]++] = facei;
288 
289  facei++;
290  }
291  }
292 
293  // Top faces
294  forAll(surfaceFaces, i)
295  {
296  eCells[i][nCellFaces[i]++] = facei;
297 
298  facei++;
299  }
300 
301  // Bottom faces
302  forAll(surfaceFaces, i)
303  {
304  label cell0 = (nLayers-1)*surfaceFaces.size() + i;
305 
306  eCells[cell0][nCellFaces[cell0]++] = facei;
307 
308  facei++;
309  }
310 
311  // return points for transferring
312  return xferMove(eCells);
313 }
314 
315 
316 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
317 
318 template
319 <
320  class Face,
321  template<class> class FaceList,
322  class PointField
323 >
324 Foam::extrudedMesh::extrudedMesh
325 (
326  const IOobject& io,
327  const PrimitivePatch<Face, FaceList, PointField>& extrudePatch,
328  const extrudeModel& model
329 )
330 :
331  polyMesh
332  (
333  io,
334  extrudedPoints(extrudePatch, model),
335  extrudedFaces(extrudePatch, model),
336  extrudedCells(extrudePatch, model)
337  ),
338  model_(model)
339 {
340  List<polyPatch*> patches(3);
341 
342  label facei = nInternalFaces();
343 
344  label sz =
345  model_.nLayers()
346  *(extrudePatch.nEdges() - extrudePatch.nInternalEdges());
347 
348  patches[0] = new wallPolyPatch
349  (
350  "sides",
351  sz,
352  facei,
353  0,
354  boundaryMesh()
355  );
356 
357  facei += sz;
358 
359  patches[1] = new polyPatch
360  (
361  "originalPatch",
362  extrudePatch.size(),
363  facei,
364  1,
365  boundaryMesh()
366  );
367 
368  facei += extrudePatch.size();
369 
370  patches[2] = new polyPatch
371  (
372  "otherSide",
373  extrudePatch.size(),
374  facei,
375  2,
376  boundaryMesh()
377  );
378 
379  addPatches(patches);
380 }
381 
382 
383 // ************************ vim: set sw=4 sts=4 et: ************************ //
384