FreeFOAM The Cross-Platform CFD Toolkit
autoPatch.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  autoPatch
26 
27 Description
28  Divides external faces into patches based on (user supplied) feature
29  angle.
30 
31 Usage
32 
33  - autoPatch [OPTIONS] <feature angle [0-180]>
34 
35  @param <feature angle [0-180]> \n
36  @todo Detailed description of argument.
37 
38  @param -overwrite \n
39  Overwrite existing data.
40 
41  @param -case <dir>\n
42  Case directory.
43 
44  @param -help \n
45  Display help message.
46 
47  @param -doc \n
48  Display Doxygen API documentation page for this application.
49 
50  @param -srcDoc \n
51  Display Doxygen source documentation page for this application.
52 
53 \*---------------------------------------------------------------------------*/
54 
55 #include <OpenFOAM/argList.H>
56 #include <OpenFOAM/polyMesh.H>
57 #include <OpenFOAM/Time.H>
61 #include <OpenFOAM/OFstream.H>
62 #include <OpenFOAM/ListOps.H>
63 
64 using namespace Foam;
65 
66 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
67 
68 // Get all feature edges.
69 void collectFeatureEdges(const boundaryMesh& bMesh, labelList& markedEdges)
70 {
71  markedEdges.setSize(bMesh.mesh().nEdges());
72 
73  label markedI = 0;
74 
75  forAll(bMesh.featureSegments(), i)
76  {
77  const labelList& segment = bMesh.featureSegments()[i];
78 
79  forAll(segment, j)
80  {
81  label featEdgeI = segment[j];
82 
83  label meshEdgeI = bMesh.featureToEdge()[featEdgeI];
84 
85  markedEdges[markedI++] = meshEdgeI;
86  }
87  }
88  markedEdges.setSize(markedI);
89 }
90 
91 
92 // Main program:
93 
94 int main(int argc, char *argv[])
95 {
97  argList::validArgs.append("feature angle[0-180]");
98  argList::validOptions.insert("overwrite", "");
99 
100 # include <OpenFOAM/setRootCase.H>
101 # include <OpenFOAM/createTime.H>
102  runTime.functionObjects().off();
103 # include <OpenFOAM/createPolyMesh.H>
104  const word oldInstance = mesh.pointsInstance();
105 
106  Info<< "Mesh read in = "
107  << runTime.cpuTimeIncrement()
108  << " s\n" << endl << endl;
109 
110 
111  //
112  // Use boundaryMesh to reuse all the featureEdge stuff in there.
113  //
114 
116 
117  scalar featureAngle(readScalar(IStringStream(args.additionalArgs()[0])()));
118  bool overwrite = args.optionFound("overwrite");
119 
120  scalar minCos = Foam::cos(featureAngle * mathematicalConstant::pi/180.0);
121 
122  Info<< "Feature:" << featureAngle << endl
123  << "minCos :" << minCos << endl
124  << endl;
125 
126  bMesh.read(mesh);
127 
128  // Set feature angle (calculate feature edges)
129  bMesh.setFeatureEdges(minCos);
130 
131  // Collect all feature edges as edge labels
132  labelList markedEdges;
133 
134  collectFeatureEdges(bMesh, markedEdges);
135 
136 
137 
138  // (new) patch ID for every face in mesh.
139  labelList patchIDs(bMesh.mesh().size(), -1);
140 
141  //
142  // Fill patchIDs with values for every face by floodfilling without
143  // crossing feature edge.
144  //
145 
146  // Current patch number.
147  label newPatchI = bMesh.patches().size();
148 
149  label suffix = 0;
150 
151  while (true)
152  {
153  // Find first unset face.
154  label unsetFaceI = findIndex(patchIDs, -1);
155 
156  if (unsetFaceI == -1)
157  {
158  // All faces have patchID set. Exit.
159  break;
160  }
161 
162  // Found unset face. Create patch for it.
163  word patchName;
164  do
165  {
166  patchName = "auto" + name(suffix++);
167  }
168  while (bMesh.findPatchID(patchName) != -1);
169 
170  bMesh.addPatch(patchName);
171 
172  bMesh.changePatchType(patchName, "patch");
173 
174 
175  // Fill visited with all faces reachable from unsetFaceI.
176  boolList visited(bMesh.mesh().size());
177 
178  bMesh.markFaces(markedEdges, unsetFaceI, visited);
179 
180 
181  // Assign all visited faces to current patch
182  label nVisited = 0;
183 
184  forAll(visited, faceI)
185  {
186  if (visited[faceI])
187  {
188  nVisited++;
189 
190  patchIDs[faceI] = newPatchI;
191  }
192  }
193 
194  Info<< "Assigned " << nVisited << " faces to patch " << patchName
195  << endl << endl;
196 
197  newPatchI++;
198  }
199 
200 
201 
202  const PtrList<boundaryPatch>& patches = bMesh.patches();
203 
204  // Create new list of patches with old ones first
205  List<polyPatch*> newPatchPtrList(patches.size());
206 
207  newPatchI = 0;
208 
209  // Copy old patches
210  forAll(mesh.boundaryMesh(), patchI)
211  {
212  const polyPatch& patch = mesh.boundaryMesh()[patchI];
213 
214  newPatchPtrList[newPatchI] =
215  patch.clone
216  (
217  mesh.boundaryMesh(),
218  newPatchI,
219  patch.size(),
220  patch.start()
221  ).ptr();
222 
223  newPatchI++;
224  }
225 
226  // Add new ones with empty size.
227  for (label patchI = newPatchI; patchI < patches.size(); patchI++)
228  {
229  const boundaryPatch& bp = patches[patchI];
230 
231  newPatchPtrList[newPatchI] = polyPatch::New
232  (
233  polyPatch::typeName,
234  bp.name(),
235  0,
236  mesh.nFaces(),
237  newPatchI,
239  ).ptr();
240 
241  newPatchI++;
242  }
243 
244  if (!overwrite)
245  {
246  runTime++;
247  }
248 
249 
250  // Change patches
251  repatchPolyTopoChanger polyMeshRepatcher(mesh);
252  polyMeshRepatcher.changePatches(newPatchPtrList);
253 
254 
255  // Change face ordering
256 
257  // Since bMesh read from mesh there is one to one mapping so we don't
258  // have to do the geometric stuff.
259  const labelList& meshFace = bMesh.meshFace();
260 
261  forAll(patchIDs, faceI)
262  {
263  label meshFaceI = meshFace[faceI];
264 
265  polyMeshRepatcher.changePatchID(meshFaceI, patchIDs[faceI]);
266  }
267 
268  polyMeshRepatcher.repatch();
269 
270  // Write resulting mesh
271  if (overwrite)
272  {
273  mesh.setInstance(oldInstance);
274  }
275  mesh.write();
276 
277 
278  Info<< "End\n" << endl;
279 
280  return 0;
281 }
282 
283 
284 // ************************ vim: set sw=4 sts=4 et: ************************ //