FreeFOAM The Cross-Platform CFD Toolkit
createBaffles.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  createBaffles
26 
27 Description
28  Makes internal faces into boundary faces.
29 
30  Does not duplicate points, unlike mergeOrSplitBaffles.
31 
32 Note
33  If any coupled patch face is selected for baffling automatically
34  the opposite member has to be selected for baffling as well. Note that this
35  is the same as repatching. This was added only for convenience so
36  you don't have to filter coupled boundary out of your set.
37 
38 Usage
39 
40  - createBaffles [OPTIONS] <cellSet name> <patchName>
41 
42  @param <cellSet name> \n
43  @todo Detailed description of argument.
44 
45  @param <patchName> \n
46  @todo Detailed description of argument.
47 
48  @param -overwrite \n
49  Overwrite existing data.
50 
51  @param -case <dir>\n
52  Case directory.
53 
54  @param -parallel \n
55  Run in parallel.
56 
57  @param -help \n
58  Display help message.
59 
60  @param -doc \n
61  Display Doxygen API documentation page for this application.
62 
63  @param -srcDoc \n
64  Display Doxygen source documentation page for this application.
65 
66 \*---------------------------------------------------------------------------*/
67 
68 #include <OpenFOAM/syncTools.H>
69 #include <OpenFOAM/argList.H>
70 #include <OpenFOAM/Time.H>
71 #include <meshTools/faceSet.H>
75 #include <OpenFOAM/ReadFields.H>
76 #include <finiteVolume/volFields.H>
78 #include <OpenFOAM/ZoneIDs.H>
79 
80 using namespace Foam;
81 
82 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
83 
84 void modifyOrAddFace
85 (
86  polyTopoChange& meshMod,
87  const face& f,
88  const label faceI,
89  const label own,
90  const bool flipFaceFlux,
91  const label newPatchI,
92  const label zoneID,
93  const bool zoneFlip,
94 
95  PackedBoolList& modifiedFace
96 )
97 {
98  if (!modifiedFace[faceI])
99  {
100  // First usage of face. Modify.
101  meshMod.setAction
102  (
104  (
105  f, // modified face
106  faceI, // label of face
107  own, // owner
108  -1, // neighbour
109  flipFaceFlux, // face flip
110  newPatchI, // patch for face
111  false, // remove from zone
112  zoneID, // zone for face
113  zoneFlip // face flip in zone
114  )
115  );
116  modifiedFace[faceI] = 1;
117  }
118  else
119  {
120  // Second or more usage of face. Add.
121  meshMod.setAction
122  (
124  (
125  f, // modified face
126  own, // owner
127  -1, // neighbour
128  -1, // master point
129  -1, // master edge
130  faceI, // master face
131  flipFaceFlux, // face flip
132  newPatchI, // patch for face
133  zoneID, // zone for face
134  zoneFlip // face flip in zone
135  )
136  );
137  }
138 }
139 
140 
141 label findPatchID(const polyMesh& mesh, const word& name)
142 {
143  label patchI = mesh.boundaryMesh().findPatchID(name);
144 
145  if (patchI == -1)
146  {
147  FatalErrorIn("findPatchID(const polyMesh&, const word&)")
148  << "Cannot find patch " << name << endl
149  << "Valid patches are " << mesh.boundaryMesh().names()
150  << exit(FatalError);
151  }
152  return patchI;
153 }
154 
155 
156 // Main program:
157 
158 int main(int argc, char *argv[])
159 {
160 # include <OpenFOAM/addRegionOption.H>
161  argList::validArgs.append("faceZone");
162  argList::validArgs.append("patch");
163  argList::validOptions.insert("additionalPatches", "(patch2 .. patchN)");
164  argList::validOptions.insert("internalFacesOnly", "");
165  argList::validOptions.insert("overwrite", "");
166 
167 # include <OpenFOAM/setRootCase.H>
168 # include <OpenFOAM/createTime.H>
169  runTime.functionObjects().off();
170 # include <OpenFOAM/createNamedMesh.H>
171  const word oldInstance = mesh.pointsInstance();
172 
173  const polyBoundaryMesh& patches = mesh.boundaryMesh();
174  const faceZoneMesh& faceZones = mesh.faceZones();
175 
176  // Faces to baffle
177  faceZoneID zoneID(args.additionalArgs()[0], faceZones);
178 
179  Info<< "Converting faces on zone " << zoneID.name()
180  << " into baffles." << nl << endl;
181 
182  if (zoneID.index() == -1)
183  {
184  FatalErrorIn(args.executable()) << "Cannot find faceZone "
185  << zoneID.name() << endl
186  << "Valid zones are " << faceZones.names()
187  << exit(FatalError);
188  }
189 
190  const faceZone& fZone = faceZones[zoneID.index()];
191 
192  Info<< "Found " << returnReduce(fZone.size(), sumOp<label>())
193  << " faces on zone " << zoneID.name() << nl << endl;
194 
195  // Make sure patches and zoneFaces are synchronised across couples
196  patches.checkParallelSync(true);
197  fZone.checkParallelSync(true);
198 
199  // Patches to put baffles into
200  DynamicList<label> newPatches(1);
201 
202  word patchName(args.additionalArgs()[1]);
203  newPatches.append(findPatchID(mesh, patchName));
204  Info<< "Using patch " << patchName
205  << " at index " << newPatches[0] << endl;
206 
207 
208  // Additional patches
209  if (args.optionFound("additionalPatches"))
210  {
211  const wordList patchNames
212  (
213  args.optionLookup("additionalPatches")()
214  );
215 
216  newPatches.reserve(patchNames.size() + 1);
217  forAll(patchNames, i)
218  {
219  newPatches.append(findPatchID(mesh, patchNames[i]));
220  Info<< "Using additional patch " << patchNames[i]
221  << " at index " << newPatches[newPatches.size()-1] << endl;
222  }
223  }
224 
225 
226  bool overwrite = args.optionFound("overwrite");
227 
228  bool internalFacesOnly = args.optionFound("internalFacesOnly");
229 
230  if (internalFacesOnly)
231  {
232  Info<< "Not converting faces on non-coupled patches." << nl << endl;
233  }
234 
235 
236  // Read objects in time directory
237  IOobjectList objects(mesh, runTime.timeName());
238 
239  // Read vol fields.
240 
242  ReadFields(mesh, objects, vsFlds);
243 
245  ReadFields(mesh, objects, vvFlds);
246 
248  ReadFields(mesh, objects, vstFlds);
249 
250  PtrList<volSymmTensorField> vsymtFlds;
251  ReadFields(mesh, objects, vsymtFlds);
252 
254  ReadFields(mesh, objects, vtFlds);
255 
256  // Read surface fields.
257 
259  ReadFields(mesh, objects, ssFlds);
260 
262  ReadFields(mesh, objects, svFlds);
263 
265  ReadFields(mesh, objects, sstFlds);
266 
268  ReadFields(mesh, objects, ssymtFlds);
269 
271  ReadFields(mesh, objects, stFlds);
272 
273 
274  // Mesh change container
275  polyTopoChange meshMod(mesh);
276 
277 
278  // Do the actual changes. Note:
279  // - loop in incrementing face order (not necessary if faceZone ordered).
280  // Preserves any existing ordering on patch faces.
281  // - two passes, do non-flip faces first and flip faces second. This
282  // guarantees that when e.g. creating a cyclic all faces from one
283  // side come first and faces from the other side next.
284 
285  // Whether first use of face (modify) or consecutive (add)
286  PackedBoolList modifiedFace(mesh.nFaces());
287  // Never modify coupled faces
288  forAll(patches, patchI)
289  {
290  const polyPatch& pp = patches[patchI];
291  if (pp.coupled())
292  {
293  forAll(pp, i)
294  {
295  modifiedFace[pp.start()+i] = 1;
296  }
297  }
298  }
299  label nModified = 0;
300 
301  forAll(newPatches, i)
302  {
303  label newPatchI = newPatches[i];
304 
305  // Pass 1. Do selected side of zone
306  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
307 
308  for (label faceI = 0; faceI < mesh.nInternalFaces(); faceI++)
309  {
310  label zoneFaceI = fZone.whichFace(faceI);
311 
312  if (zoneFaceI != -1)
313  {
314  if (!fZone.flipMap()[zoneFaceI])
315  {
316  // Use owner side of face
317  modifyOrAddFace
318  (
319  meshMod,
320  mesh.faces()[faceI], // modified face
321  faceI, // label of face
322  mesh.faceOwner()[faceI],// owner
323  false, // face flip
324  newPatchI, // patch for face
325  zoneID.index(), // zone for face
326  false, // face flip in zone
327  modifiedFace // modify or add status
328  );
329  }
330  else
331  {
332  // Use neighbour side of face
333  modifyOrAddFace
334  (
335  meshMod,
336  mesh.faces()[faceI].reverseFace(), // modified face
337  faceI, // label of face
338  mesh.faceNeighbour()[faceI],// owner
339  true, // face flip
340  newPatchI, // patch for face
341  zoneID.index(), // zone for face
342  true, // face flip in zone
343  modifiedFace // modify or add status
344  );
345  }
346 
347  nModified++;
348  }
349  }
350 
351 
352  // Pass 2. Do other side of zone
353  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
354 
355  for (label faceI = 0; faceI < mesh.nInternalFaces(); faceI++)
356  {
357  label zoneFaceI = fZone.whichFace(faceI);
358 
359  if (zoneFaceI != -1)
360  {
361  if (!fZone.flipMap()[zoneFaceI])
362  {
363  // Use neighbour side of face
364  modifyOrAddFace
365  (
366  meshMod,
367  mesh.faces()[faceI].reverseFace(), // modified face
368  faceI, // label of face
369  mesh.faceNeighbour()[faceI], // owner
370  true, // face flip
371  newPatchI, // patch for face
372  zoneID.index(), // zone for face
373  true, // face flip in zone
374  modifiedFace // modify or add
375  );
376  }
377  else
378  {
379  // Use owner side of face
380  modifyOrAddFace
381  (
382  meshMod,
383  mesh.faces()[faceI], // modified face
384  faceI, // label of face
385  mesh.faceOwner()[faceI],// owner
386  false, // face flip
387  newPatchI, // patch for face
388  zoneID.index(), // zone for face
389  false, // face flip in zone
390  modifiedFace // modify or add status
391  );
392  }
393  }
394  }
395 
396 
397  // Modify any boundary faces
398  // ~~~~~~~~~~~~~~~~~~~~~~~~~
399 
400  // Normal boundary:
401  // - move to new patch. Might already be back-to-back baffle
402  // you want to add cyclic to. Do warn though.
403  //
404  // Processor boundary:
405  // - do not move to cyclic
406  // - add normal patches though.
407 
408  // For warning once per patch.
409  labelHashSet patchWarned;
410 
411  forAll(patches, patchI)
412  {
413  const polyPatch& pp = patches[patchI];
414 
415  if (pp.coupled() && patches[newPatchI].coupled())
416  {
417  // Do not allow coupled faces to be moved to different coupled
418  // patches.
419  }
420  else if (pp.coupled() || !internalFacesOnly)
421  {
422  forAll(pp, i)
423  {
424  label faceI = pp.start()+i;
425 
426  label zoneFaceI = fZone.whichFace(faceI);
427 
428  if (zoneFaceI != -1)
429  {
430  if (patchWarned.insert(patchI))
431  {
433  << "Found boundary face (in patch " << pp.name()
434  << ") in faceZone " << fZone.name()
435  << " to convert to baffle patch "
436  << patches[newPatchI].name()
437  << endl
438  << " Run with -internalFacesOnly option"
439  << " if you don't wish to convert"
440  << " boundary faces." << endl;
441  }
442 
443  modifyOrAddFace
444  (
445  meshMod,
446  mesh.faces()[faceI], // modified face
447  faceI, // label of face
448  mesh.faceOwner()[faceI], // owner
449  false, // face flip
450  newPatchI, // patch for face
451  zoneID.index(), // zone for face
452  fZone.flipMap()[zoneFaceI], // face flip in zone
453  modifiedFace // modify or add status
454  );
455  nModified++;
456  }
457  }
458  }
459  }
460  }
461 
462 
463  Info<< "Converted " << returnReduce(nModified, sumOp<label>())
464  << " faces into boundary faces on patch " << patchName << nl << endl;
465 
466  if (!overwrite)
467  {
468  runTime++;
469  }
470 
471  // Change the mesh. Change points directly (no inflation).
472  autoPtr<mapPolyMesh> map = meshMod.changeMesh(mesh, false);
473 
474  // Update fields
475  mesh.updateMesh(map);
476 
477  // Move mesh (since morphing might not do this)
478  if (map().hasMotionPoints())
479  {
480  mesh.movePoints(map().preMotionPoints());
481  }
482 
483  if (overwrite)
484  {
485  mesh.setInstance(oldInstance);
486  }
487  Info<< "Writing mesh to " << runTime.timeName() << endl;
488 
489  mesh.write();
490 
491  Info<< "End\n" << endl;
492 
493  return 0;
494 }
495 
496 
497 // ************************ vim: set sw=4 sts=4 et: ************************ //