1 /*---------------------------------------------------------------------------*\
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
5 \\ / A nd | Copyright (C) 1991-2008 OpenCFD Ltd.
7 -------------------------------------------------------------------------------
9 This file is part of OpenFOAM.
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 the
13 Free Software Foundation; either version 2 of the License, or (at your
14 option) any later version.
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
21 You should have received a copy of the GNU General Public License
22 along with OpenFOAM; if not, write to the Free Software Foundation,
23 Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
25 \*---------------------------------------------------------------------------*/
27 #include "extrudedMesh.H"
28 #include "wallPolyPatch.H"
29 #include "meshTools.H"
33 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
35 bool Foam::extrudedMesh::sameOrder(const face& f, const edge& e)
37 label i = findIndex(f, e[0]);
39 label nextI = (i == f.size()-1 ? 0 : i+1);
41 return f[nextI] == e[1];
48 template<class> class FaceList,
51 Foam::pointField Foam::extrudedMesh::extrudedPoints
53 const PrimitivePatch<Face, FaceList, PointField>& extrudePatch,
54 const extrudeModel& model
57 const pointField& surfacePoints = extrudePatch.localPoints();
58 const vectorField& surfaceNormals = extrudePatch.pointNormals();
60 const label nLayers = model.nLayers();
62 pointField ePoints((nLayers + 1)*surfacePoints.size());
64 for (label layer=0; layer<=nLayers; layer++)
66 label offset = layer*surfacePoints.size();
68 forAll(surfacePoints, i)
70 ePoints[offset + i] = model
83 template<class Face, template<class> class FaceList, class PointField>
84 Foam::faceList Foam::extrudedMesh::extrudedFaces
86 const PrimitivePatch<Face, FaceList, PointField>& extrudePatch,
87 const extrudeModel& model
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();
95 const label nLayers = model.nLayers();
98 (nLayers + 1)*surfaceFaces.size() + nLayers*surfaceEdges.size();
100 faceList eFaces(nFaces);
106 for (label layer=0; layer<nLayers; layer++)
108 label currentLayerOffset = layer*surfacePoints.size();
109 label nextLayerOffset = currentLayerOffset + surfacePoints.size();
111 // Side faces from layer to layer+1
112 for (label i=0; i<nInternalEdges; i++)
114 quad[0] = surfaceEdges[i][1] + currentLayerOffset;
115 quad[1] = surfaceEdges[i][0] + currentLayerOffset;
116 quad[2] = surfaceEdges[i][0] + nextLayerOffset;
117 quad[3] = surfaceEdges[i][1] + nextLayerOffset;
119 eFaces[facei++] = face(quad);
122 // Faces between layer and layer+1
123 if (layer < nLayers-1)
125 forAll(surfaceFaces, i)
130 surfaceFaces[i].reverseFace()
137 // External side faces
138 for (label layer=0; layer<nLayers; layer++)
140 label currentLayerOffset = layer*surfacePoints.size();
141 label nextLayerOffset = currentLayerOffset + surfacePoints.size();
143 // Side faces across layer
144 for (label i=nInternalEdges; i<surfaceEdges.size(); i++)
146 const edge& e = surfaceEdges[i];
147 quad[0] = e[1] + currentLayerOffset;
148 quad[1] = e[0] + currentLayerOffset;
149 quad[2] = e[0] + nextLayerOffset;
150 quad[3] = e[1] + nextLayerOffset;
152 label ownerFace = extrudePatch.edgeFaces()[i][0];
154 if (!sameOrder(surfaceFaces[ownerFace], e))
159 eFaces[facei++] = face(quad);
164 forAll(surfaceFaces, i)
166 eFaces[facei++] = face(surfaceFaces[i]);
170 forAll(surfaceFaces, i)
175 surfaceFaces[i].reverseFace()
176 + nLayers*surfacePoints.size()
184 template<class Face, template<class> class FaceList, class PointField>
185 Foam::cellList Foam::extrudedMesh::extrudedCells
187 const PrimitivePatch<Face, FaceList, PointField>& extrudePatch,
188 const extrudeModel& model
191 const List<face>& surfaceFaces = extrudePatch.localFaces();
192 const edgeList& surfaceEdges = extrudePatch.edges();
193 const label nInternalEdges = extrudePatch.nInternalEdges();
195 const label nLayers = model.nLayers();
197 cellList eCells(nLayers*surfaceFaces.size());
200 forAll(surfaceFaces, i)
202 const face& f = surfaceFaces[i];
204 for (label layer=0; layer<nLayers; layer++)
206 eCells[i + layer*surfaceFaces.size()].setSize(f.size() + 2);
210 // Current face count per cell.
211 labelList nCellFaces(eCells.size(), 0);
216 for (label layer=0; layer<nLayers; layer++)
218 // Side faces from layer to layer+1
219 for (label i=0; i<nInternalEdges; i++)
221 // Get patch faces using edge
222 const labelList& edgeFaces = extrudePatch.edgeFaces()[i];
224 // Get cells on this layer
225 label cell0 = layer*surfaceFaces.size() + edgeFaces[0];
226 label cell1 = layer*surfaceFaces.size() + edgeFaces[1];
228 eCells[cell0][nCellFaces[cell0]++] = facei;
229 eCells[cell1][nCellFaces[cell1]++] = facei;
234 // Faces between layer and layer+1
235 if (layer < nLayers-1)
237 forAll(surfaceFaces, i)
239 label cell0 = layer*surfaceFaces.size() + i;
240 label cell1 = (layer+1)*surfaceFaces.size() + i;
242 eCells[cell0][nCellFaces[cell0]++] = facei;
243 eCells[cell1][nCellFaces[cell1]++] = facei;
250 // External side faces
251 for (label layer=0; layer<nLayers; layer++)
253 // Side faces across layer
254 for (label i=nInternalEdges; i<surfaceEdges.size(); i++)
256 // Get patch faces using edge
257 const labelList& edgeFaces = extrudePatch.edgeFaces()[i];
259 // Get cells on this layer
260 label cell0 = layer*surfaceFaces.size() + edgeFaces[0];
262 eCells[cell0][nCellFaces[cell0]++] = facei;
269 forAll(surfaceFaces, i)
271 eCells[i][nCellFaces[i]++] = facei;
277 forAll(surfaceFaces, i)
279 label cell0 = (nLayers-1)*surfaceFaces.size() + i;
281 eCells[cell0][nCellFaces[cell0]++] = facei;
290 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
295 template<class> class FaceList,
298 Foam::extrudedMesh::extrudedMesh
301 const PrimitivePatch<Face, FaceList, PointField>& extrudePatch,
302 const extrudeModel& model
308 extrudedPoints(extrudePatch, model),
309 extrudedFaces(extrudePatch, model),
310 extrudedCells(extrudePatch, model)
314 List<polyPatch*> patches(3);
316 label facei = nInternalFaces();
320 *(extrudePatch.nEdges() - extrudePatch.nInternalEdges());
322 patches[0] = new wallPolyPatch
333 patches[1] = new polyPatch
342 facei += extrudePatch.size();
344 patches[2] = new polyPatch
357 // ************************************************************************* //