FreeFOAM The Cross-Platform CFD Toolkit
surfaceRedistributePar.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-2007 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  surfaceRedistributePar
26 
27 Description
28  (Re)distribution of triSurface. Either takes an undecomposed surface
29  or an already decomposed surface and redistribute it so each processor
30  has all triangles that overlap its mesh.
31 
32 Usage
33  - surfaceRedistributePar [OPTION]
34 
35  @param <triSurfaceMesh> \n
36  Tri-Surface mesh to read.
37 
38  @param <distributionType> \n
39  @todo Detailed description.
40 
41  @param -keepNonMapped \n
42  Preserve surface outside of mesh bounds.
43 
44  @param -case <dir> \n
45  Specify the case directory
46 
47  @param -parallel \n
48  Run the case in parallel
49 
50  @param -help \n
51  Display short usage message
52 
53  @param -doc \n
54  Display Doxygen documentation page
55 
56  @param -srcDoc \n
57  Display source code
58 
59 Note
60  - best decomposition option is hierarchGeomDecomp since
61  guarantees square decompositions.
62  - triangles might be present on multiple processors.
63  - merging uses geometric tolerance so take care with writing precision.
64 
65 \*---------------------------------------------------------------------------*/
66 
67 #include <meshTools/treeBoundBox.H>
68 #include <OpenFOAM/FixedList.H>
69 #include <OpenFOAM/argList.H>
70 #include <OpenFOAM/Time.H>
71 #include <OpenFOAM/polyMesh.H>
73 #include <OpenFOAM/mapDistribute.H>
75 #include <OpenFOAM/Pair.H>
76 
77 using namespace Foam;
78 
79 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
80 
81 // Print on master all the per-processor surface stats.
82 void writeProcStats
83 (
84  const triSurface& s,
85  const List<List<treeBoundBox> >& meshBb
86 )
87 {
88  // Determine surface bounding boxes, faces, points
90  {
91  surfBb[Pstream::myProcNo()] = treeBoundBox(s.points());
92  Pstream::gatherList(surfBb);
93  Pstream::scatterList(surfBb);
94  }
95 
100 
101  labelList nFaces(Pstream::nProcs());
102  nFaces[Pstream::myProcNo()] = s.size();
103  Pstream::gatherList(nFaces);
104  Pstream::scatterList(nFaces);
105 
106  forAll(surfBb, procI)
107  {
108  const List<treeBoundBox>& bbs = meshBb[procI];
109 
110  Info<< "processor" << procI << endl
111  << "\tMesh bounds : " << bbs[0] << nl;
112  for (label i = 1; i < bbs.size(); i++)
113  {
114  Info<< "\t " << bbs[i]<< nl;
115  }
116  Info<< "\tSurface bounding box : " << surfBb[procI] << nl
117  << "\tTriangles : " << nFaces[procI] << nl
118  << "\tVertices : " << nPoints[procI]
119  << endl;
120  }
121  Info<< endl;
122 }
123 
124 
125 // Main program:
126 
127 int main(int argc, char *argv[])
128 {
129  argList::validArgs.append("triSurfaceMesh");
130  argList::validArgs.append("distributionType");
131 
132  argList::validOptions.insert("keepNonMapped", "");
133 # include <OpenFOAM/setRootCase.H>
134 # include <OpenFOAM/createTime.H>
135  runTime.functionObjects().off();
136 
137  fileName surfFileName(args.additionalArgs()[0]);
138  Info<< "Reading surface from " << surfFileName << nl << endl;
139 
140  const word distType(args.additionalArgs()[1]);
141 
142  Info<< "Using distribution method "
144  << " " << distType << nl << endl;
145 
146  bool keepNonMapped = args.options().found("keepNonMapped");
147 
148  if (keepNonMapped)
149  {
150  Info<< "Preserving surface outside of mesh bounds." << nl << endl;
151  }
152  else
153  {
154  Info<< "Removing surface outside of mesh bounds." << nl << endl;
155  }
156 
157 
158  if (!Pstream::parRun())
159  {
161  << "Please run this program on the decomposed case."
162  << " It will read surface " << surfFileName
163  << " and decompose it such that it overlaps the mesh bounding box."
164  << exit(FatalError);
165  }
166 
167 
168 # include <OpenFOAM/createPolyMesh.H>
169 
170  Random rndGen(653213);
171 
172  // Determine mesh bounding boxes:
174  {
176  (
177  1,
179  (
180  boundBox(mesh.points(), false)
181  ).extend(rndGen, 1E-3)
182  );
183  Pstream::gatherList(meshBb);
184  Pstream::scatterList(meshBb);
185  }
186 
187  IOobject io
188  (
189  surfFileName, // name
190  //runTime.findInstance("triSurface", surfFileName), // instance
191  runTime.constant(), // instance
192  "triSurface", // local
193  runTime, // registry
196  );
197 
198  const fileName actualPath(io.filePath());
199  fileName localPath(actualPath);
200  localPath.replace(runTime.rootPath() + '/', "");
201 
202  if (actualPath == io.objectPath())
203  {
204  Info<< "Loading local (decomposed) surface " << localPath << nl <<endl;
205  }
206  else
207  {
208  Info<< "Loading undecomposed surface " << localPath << nl << endl;
209  }
210 
211 
212  // Create dummy dictionary for bounding boxes if does not exist.
213  if (!isFile(actualPath / "Dict"))
214  {
215  dictionary dict;
216  dict.add("bounds", meshBb[Pstream::myProcNo()]);
217  dict.add("distributionType", distType);
218  dict.add("mergeDistance", SMALL);
219 
220  IOdictionary ioDict
221  (
222  IOobject
223  (
224  io.name() + "Dict",
225  io.instance(),
226  io.local(),
227  io.db(),
230  false
231  ),
232  dict
233  );
234 
235  Info<< "Writing dummy bounds dictionary to " << ioDict.name()
236  << nl << endl;
237 
238  ioDict.regIOobject::writeObject
239  (
242  ioDict.time().writeCompression()
243  );
244  }
245 
246 
247  // Load surface
249  Info<< "Loaded surface" << nl << endl;
250 
251 
252  // Generate a test field
253  {
254  const triSurface& s = static_cast<const triSurface&>(surfMesh);
255 
257  (
259  (
260  IOobject
261  (
262  surfMesh.searchableSurface::name(), // name
263  surfMesh.searchableSurface::instance(), // instance
264  surfMesh.searchableSurface::local(), // local
265  surfMesh,
268  ),
269  surfMesh,
270  dimLength
271  )
272  );
273  triSurfaceVectorField& fc = fcPtr();
274 
275  forAll(fc, triI)
276  {
277  fc[triI] = s[triI].centre(s.points());
278  }
279 
280  // Steal pointer and store object on surfMesh
281  fcPtr.ptr()->store();
282  }
283 
284 
285  // Write per-processor stats
286  Info<< "Before redistribution:" << endl;
287  writeProcStats(surfMesh, meshBb);
288 
289 
290  // Do redistribution
291  Info<< "Redistributing surface" << nl << endl;
292  autoPtr<mapDistribute> faceMap;
293  autoPtr<mapDistribute> pointMap;
294  surfMesh.distribute
295  (
296  meshBb[Pstream::myProcNo()],
297  keepNonMapped,
298  faceMap,
299  pointMap
300  );
301  faceMap.clear();
302  pointMap.clear();
303 
304  Info<< endl;
305 
306 
307  // Write per-processor stats
308  Info<< "After redistribution:" << endl;
309  writeProcStats(surfMesh, meshBb);
310 
311 
312  Info<< "Writing surface." << nl << endl;
313  surfMesh.searchableSurface::write();
314 
315  Info<< "End\n" << endl;
316 
317  return 0;
318 }
319 
320 
321 // ************************ vim: set sw=4 sts=4 et: ************************ //