FreeFOAM The Cross-Platform CFD Toolkit
extrudeMesh.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 Application
25  extrudeMesh
26 
27 Description
28  Extrude mesh from existing patch or from patch read from file.
29 
30  By default uses outwards facing normals, optionally flips faces.
31 
32  Type of extrusion prescribed by run-time selectable model.
33 
34 Usage
35 
36  - extrudeMesh [OPTIONS]
37 
38  @param -mergeFaces \n
39  Merge faces for axisymmetric cases.
40 
41  @param -sourceRoot <directory containing the case>\n
42  Parent directory of the source case.
43 
44  @param -sourceCase <source case name>\n
45  Source case name.
46 
47  @param -sourcePatch <name of source patch to extrude>\n
48  Source patch.
49 
50  @param -surface <Foam surface file>\n
51  Specify surface file to extrude.
52 
53  @param -case <dir>\n
54  Case directory.
55 
56  @param -help \n
57  Display help message.
58 
59  @param -doc \n
60  Display Doxygen API documentation page for this application.
61 
62  @param -srcDoc \n
63  Display Doxygen source documentation page for this application.
64 
65 Note
66  Merges close points so be careful.
67 
68 \*---------------------------------------------------------------------------*/
69 
70 #include <OpenFOAM/argList.H>
71 #include <OpenFOAM/Time.H>
73 #include <OpenFOAM/IFstream.H>
74 #include "faceMesh.H"
81 
84 
85 using namespace Foam;
86 
87 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
88 // Main program:
89 
90 int main(int argc, char *argv[])
91 {
92  #include <OpenFOAM/setRootCase.H>
93  #include "createTimeExtruded.H"
94 
96 
97  IOdictionary dict
98  (
99  IOobject
100  (
101  "extrudeProperties",
102  runTimeExtruded.constant(),
103  runTimeExtruded,
105  )
106  );
107 
109 
110  const word sourceType(dict.lookup("constructFrom"));
111 
112  autoPtr<faceMesh> fMesh;
113 
114  if (sourceType == "patch")
115  {
116  fileName sourceCasePath(dict.lookup("sourceCase"));
117  sourceCasePath.expand();
118  fileName sourceRootDir = sourceCasePath.path();
119  fileName sourceCaseDir = sourceCasePath.name();
120  word patchName(dict.lookup("sourcePatch"));
121 
122  Info<< "Extruding patch " << patchName
123  << " on mesh " << sourceCasePath << nl
124  << endl;
125 
126  Time runTime
127  (
129  sourceRootDir,
130  sourceCaseDir
131  );
132  #include <OpenFOAM/createPolyMesh.H>
133 
134  label patchID = mesh.boundaryMesh().findPatchID(patchName);
135 
136  if (patchID == -1)
137  {
139  << "Cannot find patch " << patchName
140  << " in the source mesh.\n"
141  << "Valid patch names are " << mesh.boundaryMesh().names()
142  << exit(FatalError);
143  }
144 
145  const polyPatch& pp = mesh.boundaryMesh()[patchID];
146  fMesh.reset(new faceMesh(pp.localFaces(), pp.localPoints()));
147 
148  {
149  fileName surfName(runTime.path()/patchName + ".sMesh");
150  Info<< "Writing patch as surfaceMesh to "
151  << surfName << nl << endl;
152  OFstream os(surfName);
153  os << fMesh() << nl;
154  }
155  }
156  else if (sourceType == "surface")
157  {
158  // Read from surface
159  fileName surfName(dict.lookup("surface"));
160 
161  Info<< "Extruding surfaceMesh read from file " << surfName << nl
162  << endl;
163 
164  IFstream is(surfName);
165 
166  fMesh.reset(new faceMesh(is));
167 
168  Info<< "Read patch from file " << surfName << nl
169  << endl;
170  }
171  else
172  {
174  << "Illegal 'constructFrom' specification. Should either be "
175  << "patch or surface." << exit(FatalError);
176  }
177 
178  Switch flipNormals(dict.lookup("flipNormals"));
179 
180  if (flipNormals)
181  {
182  Info<< "Flipping faces." << nl << endl;
183 
184  faceList faces(fMesh().size());
185  forAll(faces, i)
186  {
187  faces[i] = fMesh()[i].reverseFace();
188  }
189  fMesh.reset(new faceMesh(faces, fMesh().localPoints()));
190  }
191 
192 
193  Info<< "Extruding patch with :" << nl
194  << " points : " << fMesh().points().size() << nl
195  << " faces : " << fMesh().size() << nl
196  << " normals[0] : " << fMesh().faceNormals()[0]
197  << nl
198  << endl;
199 
201  (
202  IOobject
203  (
205  runTimeExtruded.constant(),
206  runTimeExtruded
207  ),
208  fMesh(),
209  model()
210  );
211 
212 
213  const boundBox& bb = mesh.globalData().bb();
214  const vector span = bb.span();
215  const scalar mergeDim = 1E-4 * bb.minDim();
216 
217  Info<< "Mesh bounding box : " << bb << nl
218  << " with span : " << span << nl
219  << "Merge distance : " << mergeDim << nl
220  << endl;
221 
223 
224  const label origPatchID = patches.findPatchID("originalPatch");
225  const label otherPatchID = patches.findPatchID("otherSide");
226 
227  if (origPatchID == -1 || otherPatchID == -1)
228  {
230  << "Cannot find patch originalPatch or otherSide." << nl
231  << "Valid patches are " << patches.names() << exit(FatalError);
232  }
233 
234  // Collapse edges
235  // ~~~~~~~~~~~~~~
236 
237  {
238  Info<< "Collapsing edges < " << mergeDim << " ..." << nl << endl;
239 
240  // Edge collapsing engine
241  edgeCollapser collapser(mesh);
242 
243  const edgeList& edges = mesh.edges();
244  const pointField& points = mesh.points();
245 
246  forAll(edges, edgeI)
247  {
248  const edge& e = edges[edgeI];
249 
250  scalar d = e.mag(points);
251 
252  if (d < mergeDim)
253  {
254  Info<< "Merging edge " << e << " since length " << d
255  << " << " << mergeDim << nl;
256 
257  // Collapse edge to e[0]
258  collapser.collapseEdge(edgeI, e[0]);
259  }
260  }
261 
262  // Topo change container
263  polyTopoChange meshMod(mesh);
264  // Put all modifications into meshMod
265  bool anyChange = collapser.setRefinement(meshMod);
266 
267  if (anyChange)
268  {
269  // Construct new mesh from polyTopoChange.
270  autoPtr<mapPolyMesh> map = meshMod.changeMesh(mesh, false);
271 
272  // Update fields
273  mesh.updateMesh(map);
274 
275  // Move mesh (if inflation used)
276  if (map().hasMotionPoints())
277  {
278  mesh.movePoints(map().preMotionPoints());
279  }
280  }
281  }
282 
283 
284  // Merging front and back patch faces
285  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
286 
287  Switch mergeFaces(dict.lookup("mergeFaces"));
288  if (mergeFaces)
289  {
290  Info<< "Assuming full 360 degree axisymmetric case;"
291  << " stitching faces on patches "
292  << patches[origPatchID].name() << " and "
293  << patches[otherPatchID].name() << " together ..." << nl << endl;
294 
295  polyTopoChanger stitcher(mesh);
296  stitcher.setSize(1);
297 
298  // Make list of masterPatch faces
299  labelList isf(patches[origPatchID].size());
300 
301  forAll (isf, i)
302  {
303  isf[i] = patches[origPatchID].start() + i;
304  }
305 
306  const word cutZoneName("originalCutFaceZone");
307 
308  List<faceZone*> fz
309  (
310  1,
311  new faceZone
312  (
313  cutZoneName,
314  isf,
315  boolList(isf.size(), false),
316  0,
317  mesh.faceZones()
318  )
319  );
320 
322 
323  // Add the perfect interface mesh modifier
324  stitcher.set
325  (
326  0,
327  new perfectInterface
328  (
329  "couple",
330  0,
331  stitcher,
332  cutZoneName,
333  patches[origPatchID].name(),
334  patches[otherPatchID].name()
335  )
336  );
337 
338  // Execute all polyMeshModifiers
339  autoPtr<mapPolyMesh> morphMap = stitcher.changeMesh(true);
340 
341  mesh.movePoints(morphMap->preMotionPoints());
342  }
343 
344  if (!mesh.write())
345  {
346  FatalErrorIn(args.executable()) << "Failed writing mesh"
347  << exit(FatalError);
348  }
349 
350  Info << "End\n" << endl;
351 
352  return 0;
353 }
354 
355 
356 // ************************ vim: set sw=4 sts=4 et: ************************ //