FreeFOAM The Cross-Platform CFD Toolkit
mergeOrSplitBaffles.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  mergeOrSplitBaffles
26 
27 Description
28  Detects faces that share points (baffles). Either merge them or
29  duplicate the points.
30 
31 Note
32  - can only handle pairwise boundary faces. So three faces using
33  the same points is not handled (is illegal mesh anyway)
34 
35  - there is no option to only split/merge some baffles.
36 
37  - surfaces consisting of duplicate faces can be topologically split
38  if the points on the interior of the surface cannot walk to all the
39  cells that use them in one go.
40 
41  - Parallel operation (where duplicate face is perpendicular to a coupled
42  boundary) is supported but not really tested.
43  (Note that coupled faces themselves are not seen as duplicate faces)
44 
45 Usage
46 
47  - mergeOrSplitBaffles [OPTIONS]
48 
49  @param -split \n
50  Split duplicate surfaces.
51 
52  @param -overwrite \n
53  Overwrite existing data.
54 
55  @param -detectOnly \n
56  Do no processing.
57 
58  @param -case <dir>\n
59  Case directory.
60 
61  @param -parallel \n
62  Run in parallel.
63 
64  @param -help \n
65  Display help message.
66 
67  @param -doc \n
68  Display Doxygen API documentation page for this application.
69 
70  @param -srcDoc \n
71  Display Doxygen source documentation page for this application.
72 
73 \*---------------------------------------------------------------------------*/
74 
75 #include <OpenFOAM/argList.H>
76 #include <OpenFOAM/Time.H>
77 #include <OpenFOAM/syncTools.H>
78 #include <meshTools/faceSet.H>
79 #include <meshTools/pointSet.H>
80 #include <meshTools/meshTools.H>
88 #include <OpenFOAM/ReadFields.H>
89 #include <finiteVolume/volFields.H>
91 
92 using namespace Foam;
93 
94 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
95 
96 void insertDuplicateMerge
97 (
98  const polyMesh& mesh,
99  const labelList& duplicates,
100  polyTopoChange& meshMod
101 )
102 {
103  const faceList& faces = mesh.faces();
104  const labelList& faceOwner = mesh.faceOwner();
105  const faceZoneMesh& faceZones = mesh.faceZones();
106 
107  forAll(duplicates, bFaceI)
108  {
109  label otherFaceI = duplicates[bFaceI];
110 
111  if (otherFaceI != -1 && otherFaceI > bFaceI)
112  {
113  // Two duplicate faces. Merge.
114 
115  label face0 = mesh.nInternalFaces() + bFaceI;
116  label face1 = mesh.nInternalFaces() + otherFaceI;
117 
118  label own0 = faceOwner[face0];
119  label own1 = faceOwner[face1];
120 
121  if (own0 < own1)
122  {
123  // Use face0 as the new internal face.
124  label zoneID = faceZones.whichZone(face0);
125  bool zoneFlip = false;
126 
127  if (zoneID >= 0)
128  {
129  const faceZone& fZone = faceZones[zoneID];
130  zoneFlip = fZone.flipMap()[fZone.whichFace(face0)];
131  }
132 
133  meshMod.setAction(polyRemoveFace(face1));
134  meshMod.setAction
135  (
137  (
138  faces[face0], // modified face
139  face0, // label of face being modified
140  own0, // owner
141  own1, // neighbour
142  false, // face flip
143  -1, // patch for face
144  false, // remove from zone
145  zoneID, // zone for face
146  zoneFlip // face flip in zone
147  )
148  );
149  }
150  else
151  {
152  // Use face1 as the new internal face.
153  label zoneID = faceZones.whichZone(face1);
154  bool zoneFlip = false;
155 
156  if (zoneID >= 0)
157  {
158  const faceZone& fZone = faceZones[zoneID];
159  zoneFlip = fZone.flipMap()[fZone.whichFace(face1)];
160  }
161 
162  meshMod.setAction(polyRemoveFace(face0));
163  meshMod.setAction
164  (
166  (
167  faces[face1], // modified face
168  face1, // label of face being modified
169  own1, // owner
170  own0, // neighbour
171  false, // face flip
172  -1, // patch for face
173  false, // remove from zone
174  zoneID, // zone for face
175  zoneFlip // face flip in zone
176  )
177  );
178  }
179  }
180  }
181 }
182 
183 
184 labelList findBaffles(const polyMesh& mesh, const labelList& boundaryFaces)
185 {
186  // Get all duplicate face labels (in boundaryFaces indices!).
188  (
189  mesh,
190  boundaryFaces
191  );
192 
193 
194  // Check that none are on processor patches
195  const polyBoundaryMesh& patches = mesh.boundaryMesh();
196 
197  forAll(duplicates, bFaceI)
198  {
199  if (duplicates[bFaceI] != -1)
200  {
201  label faceI = mesh.nInternalFaces() + bFaceI;
202  label patchI = patches.whichPatch(faceI);
203 
204  if (isA<processorPolyPatch>(patches[patchI]))
205  {
206  FatalErrorIn("findBaffles(const polyMesh&, const labelList&)")
207  << "Duplicate face " << faceI
208  << " is on a processorPolyPatch."
209  << "This is not allowed." << nl
210  << "Face:" << faceI
211  << " is on patch:" << patches[patchI].name()
212  << abort(FatalError);
213  }
214  }
215  }
216 
217 
218  // Write to faceSet for ease of postprocessing.
219  {
220  faceSet duplicateSet
221  (
222  mesh,
223  "duplicateFaces",
224  (mesh.nFaces() - mesh.nInternalFaces())/256
225  );
226 
227  forAll(duplicates, bFaceI)
228  {
229  label otherFaceI = duplicates[bFaceI];
230 
231  if (otherFaceI != -1 && otherFaceI > bFaceI)
232  {
233  duplicateSet.insert(mesh.nInternalFaces() + bFaceI);
234  duplicateSet.insert(mesh.nInternalFaces() + otherFaceI);
235  }
236  }
237 
238  Pout<< "Writing " << duplicateSet.size()
239  << " duplicate faces to faceSet " << duplicateSet.objectPath()
240  << nl << endl;
241  duplicateSet.write();
242  }
243 
244  return duplicates;
245 }
246 
247 
248 
249 
250 int main(int argc, char *argv[])
251 {
252 # include <OpenFOAM/addRegionOption.H>
253  argList::validOptions.insert("split", "");
254  argList::validOptions.insert("overwrite", "");
255  argList::validOptions.insert("detectOnly", "");
256 # include <OpenFOAM/setRootCase.H>
257 # include <OpenFOAM/createTime.H>
258  runTime.functionObjects().off();
259 # include <OpenFOAM/createNamedMesh.H>
260  const word oldInstance = mesh.pointsInstance();
261 
262  bool split = args.optionFound("split");
263  bool overwrite = args.optionFound("overwrite");
264  bool detectOnly = args.optionFound("detectOnly");
265 
266  // Collect all boundary faces
267  labelList boundaryFaces(mesh.nFaces() - mesh.nInternalFaces());
268 
269  forAll(boundaryFaces, i)
270  {
271  boundaryFaces[i] = i+mesh.nInternalFaces();
272  }
273 
274 
275  if (detectOnly)
276  {
277  findBaffles(mesh, boundaryFaces);
278 
279  return 0;
280  }
281 
282 
283 
284  // Read objects in time directory
285  IOobjectList objects(mesh, runTime.timeName());
286 
287  // Read vol fields.
288 
290  ReadFields(mesh, objects, vsFlds);
291 
293  ReadFields(mesh, objects, vvFlds);
294 
296  ReadFields(mesh, objects, vstFlds);
297 
298  PtrList<volSymmTensorField> vsymtFlds;
299  ReadFields(mesh, objects, vsymtFlds);
300 
302  ReadFields(mesh, objects, vtFlds);
303 
304  // Read surface fields.
305 
307  ReadFields(mesh, objects, ssFlds);
308 
310  ReadFields(mesh, objects, svFlds);
311 
313  ReadFields(mesh, objects, sstFlds);
314 
316  ReadFields(mesh, objects, ssymtFlds);
317 
319  ReadFields(mesh, objects, stFlds);
320 
321 
322  // Mesh change engine
323  polyTopoChange meshMod(mesh);
324 
325 
326  if (split)
327  {
328  Pout<< "Topologically splitting duplicate surfaces"
329  << ", i.e. duplicating points internal to duplicate surfaces."
330  << nl << endl;
331 
332  // Analyse which points need to be duplicated
334 
335  // Point duplication engine
336  duplicatePoints pointDuplicator(mesh);
337 
338  // Insert topo changes
339  pointDuplicator.setRefinement(regionSide, meshMod);
340  }
341  else
342  {
343  Pout<< "Merging duplicate faces."
344  << nl << endl;
345 
346  // Get all duplicate face labels (in boundaryFaces indices!).
347  labelList duplicates(findBaffles(mesh, boundaryFaces));
348 
349  // Merge into internal faces.
350  insertDuplicateMerge(mesh, duplicates, meshMod);
351  }
352 
353  if (!overwrite)
354  {
355  runTime++;
356  }
357 
358  // Change the mesh. No inflation.
359  autoPtr<mapPolyMesh> map = meshMod.changeMesh(mesh, false);
360 
361  // Update fields
362  mesh.updateMesh(map);
363 
364  // Move mesh (since morphing does not do this)
365  if (map().hasMotionPoints())
366  {
367  mesh.movePoints(map().preMotionPoints());
368  }
369 
370  if (overwrite)
371  {
372  mesh.setInstance(oldInstance);
373  }
374  Pout<< "Writing mesh to time " << runTime.timeName() << endl;
375  mesh.write();
376 
377  // Dump duplicated points (if any)
378  if (split)
379  {
380  const labelList& pointMap = map().pointMap();
381 
382  labelList nDupPerPoint(map().nOldPoints(), 0);
383 
384  pointSet dupPoints(mesh, "duplicatedPoints", 100);
385 
386  forAll(pointMap, pointI)
387  {
388  label oldPointI = pointMap[pointI];
389 
390  nDupPerPoint[oldPointI]++;
391 
392  if (nDupPerPoint[oldPointI] > 1)
393  {
394  dupPoints.insert(map().reversePointMap()[oldPointI]);
395  dupPoints.insert(pointI);
396  }
397  }
398 
399  Pout<< "Writing " << dupPoints.size()
400  << " duplicated points to pointSet "
401  << dupPoints.objectPath() << nl << endl;
402 
403  dupPoints.write();
404  }
405 
406  Info<< "End\n" << endl;
407 
408  return 0;
409 }
410 
411 
412 // ************************ vim: set sw=4 sts=4 et: ************************ //