Merge branch 'master' of ssh://git.code.sf.net/p/foam-extend/foam-extend-3.2
[foam-extend-3.2.git] / applications / utilities / surface / surfaceRedistributePar / surfaceRedistributePar.C
blobee6fc9bd43a371c84198c9923a18c5cf40e702dc
1 /*---------------------------------------------------------------------------*\
2   =========                 |
3   \\      /  F ield         | foam-extend: Open Source CFD
4    \\    /   O peration     | Version:     3.2
5     \\  /    A nd           | Web:         http://www.foam-extend.org
6      \\/     M anipulation  | For copyright notice see file Copyright
7 -------------------------------------------------------------------------------
8 License
9     This file is part of foam-extend.
11     foam-extend is free software: you can redistribute it and/or modify it
12     under the terms of the GNU General Public License as published by the
13     Free Software Foundation, either version 3 of the License, or (at your
14     option) any later version.
16     foam-extend is distributed in the hope that it will be useful, but
17     WITHOUT ANY WARRANTY; without even the implied warranty of
18     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
19     General Public License for more details.
21     You should have received a copy of the GNU General Public License
22     along with foam-extend.  If not, see <http://www.gnu.org/licenses/>.
24 Application
25     surfaceRedistributePar
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.
32 Note
33     - best decomposition option is hierarchGeomDecomp since
34     guarantees square decompositions.
35     - triangles might be present on multiple processors.
36     - merging uses geometric tolerance so take care with writing precision.
38 \*---------------------------------------------------------------------------*/
40 #include "treeBoundBox.H"
41 #include "FixedList.H"
42 #include "argList.H"
43 #include "objectRegistry.H"
44 #include "foamTime.H"
45 #include "polyMesh.H"
46 #include "distributedTriSurfaceMesh.H"
47 #include "mapDistribute.H"
48 #include "triSurfaceFields.H"
49 #include "Pair.H"
51 using namespace Foam;
53 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
55 // Print on master all the per-processor surface stats.
56 void writeProcStats
58     const triSurface& s,
59     const List<List<treeBoundBox> >& meshBb
62     // Determine surface bounding boxes, faces, points
63     List<treeBoundBox> surfBb(Pstream::nProcs());
64     {
65         surfBb[Pstream::myProcNo()] = treeBoundBox(s.points());
66         Pstream::gatherList(surfBb);
67         Pstream::scatterList(surfBb);
68     }
70     labelList nPoints(Pstream::nProcs());
71     nPoints[Pstream::myProcNo()] = s.points().size();
72     Pstream::gatherList(nPoints);
73     Pstream::scatterList(nPoints);
75     labelList nFaces(Pstream::nProcs());
76     nFaces[Pstream::myProcNo()] = s.size();
77     Pstream::gatherList(nFaces);
78     Pstream::scatterList(nFaces);
80     forAll(surfBb, procI)
81     {
82         const List<treeBoundBox>& bbs = meshBb[procI];
84         Info<< "processor" << procI << endl
85             << "\tMesh bounds          : " << bbs[0] << nl;
86         for (label i = 1; i < bbs.size(); i++)
87         {
88             Info<< "\t                       " << bbs[i]<< nl;
89         }
90         Info<< "\tSurface bounding box : " << surfBb[procI] << nl
91             << "\tTriangles            : " << nFaces[procI] << nl
92             << "\tVertices             : " << nPoints[procI]
93             << endl;
94     }
95     Info<< endl;
99 // Main program:
101 int main(int argc, char *argv[])
103     argList::validArgs.append("triSurfaceMesh");
104     argList::validArgs.append("distributionType");
106     argList::validOptions.insert("keepNonMapped", "");
107 #   include "setRootCase.H"
108 #   include "createTime.H"
109     runTime.functionObjects().off();
111     fileName surfFileName(args.additionalArgs()[0]);
112     Info<< "Reading surface from " << surfFileName << nl << endl;
114     const word distType(args.additionalArgs()[1]);
116     Info<< "Using distribution method "
117         << distributedTriSurfaceMesh::distributionTypeNames_[distType]
118         << " " << distType << nl << endl;
120     bool keepNonMapped = args.options().found("keepNonMapped");
122     if (keepNonMapped)
123     {
124         Info<< "Preserving surface outside of mesh bounds." << nl << endl;
125     }
126     else
127     {
128         Info<< "Removing surface outside of mesh bounds." << nl << endl;
129     }
132     if (!Pstream::parRun())
133     {
134         FatalErrorIn(args.executable())
135             << "Please run this program on the decomposed case."
136             << " It will read surface " << surfFileName
137             << " and decompose it such that it overlaps the mesh bounding box."
138             << exit(FatalError);
139     }
142 #   include "createPolyMesh.H"
144     Random rndGen(653213);
146     // Determine mesh bounding boxes:
147     List<List<treeBoundBox> > meshBb(Pstream::nProcs());
148     {
149         meshBb[Pstream::myProcNo()] = List<treeBoundBox>
150         (
151             1,
152             treeBoundBox
153             (
154                 boundBox(mesh.points(), false)
155             ).extend(rndGen, 1E-3)
156         );
157         Pstream::gatherList(meshBb);
158         Pstream::scatterList(meshBb);
159     }
161     IOobject io
162     (
163         surfFileName,         // name
164         //runTime.findInstance("triSurface", surfFileName),   // instance
165         runTime.constant(),   // instance
166         "triSurface",         // local
167         runTime,              // registry
168         IOobject::MUST_READ,
169         IOobject::NO_WRITE
170     );
172     const fileName actualPath(io.filePath());
173     fileName localPath(actualPath);
174     localPath.replace(runTime.rootPath() + '/', "");
176     if (actualPath == io.objectPath())
177     {
178         Info<< "Loading local (decomposed) surface " << localPath << nl <<endl;
179     }
180     else
181     {
182         Info<< "Loading undecomposed surface " << localPath << nl << endl;
183     }
186     // Create dummy dictionary for bounding boxes if does not exist.
187     if (!isFile(actualPath / "Dict"))
188     {
189         dictionary dict;
190         dict.add("bounds", meshBb[Pstream::myProcNo()]);
191         dict.add("distributionType", distType);
192         dict.add("mergeDistance", SMALL);
194         IOdictionary ioDict
195         (
196             IOobject
197             (
198                 io.name() + "Dict",
199                 io.instance(),
200                 io.local(),
201                 io.db(),
202                 IOobject::NO_READ,
203                 IOobject::NO_WRITE,
204                 false
205             ),
206             dict
207         );
209         Info<< "Writing dummy bounds dictionary to " << ioDict.name()
210             << nl << endl;
212         ioDict.regIOobject::writeObject
213         (
214             IOstream::ASCII,
215             IOstream::currentVersion,
216             ioDict.time().writeCompression()
217         );
218     }
221     // Load surface
222     distributedTriSurfaceMesh surfMesh(io);
223     Info<< "Loaded surface" << nl << endl;
226     // Generate a test field
227     {
228         const triSurface& s = static_cast<const triSurface&>(surfMesh);
230         autoPtr<triSurfaceVectorField> fcPtr
231         (
232             new triSurfaceVectorField
233             (
234                 IOobject
235                 (
236                     surfMesh.searchableSurface::name(),     // name
237                     surfMesh.searchableSurface::instance(), // instance
238                     surfMesh.searchableSurface::local(),    // local
239                     surfMesh,
240                     IOobject::NO_READ,
241                     IOobject::AUTO_WRITE
242                 ),
243                 surfMesh,
244                 dimLength
245             )
246         );
247         triSurfaceVectorField& fc = fcPtr();
249         forAll(fc, triI)
250         {
251             fc[triI] = s[triI].centre(s.points());
252         }
254         // Steal pointer and store object on surfMesh
255         fcPtr.ptr()->store();
256     }
259     // Write per-processor stats
260     Info<< "Before redistribution:" << endl;
261     writeProcStats(surfMesh, meshBb);
264     // Do redistribution
265     Info<< "Redistributing surface" << nl << endl;
266     autoPtr<mapDistribute> faceMap;
267     autoPtr<mapDistribute> pointMap;
268     surfMesh.distribute
269     (
270         meshBb[Pstream::myProcNo()],
271         keepNonMapped,
272         faceMap,
273         pointMap
274     );
275     faceMap.clear();
276     pointMap.clear();
278     Info<< endl;
281     // Write per-processor stats
282     Info<< "After redistribution:" << endl;
283     writeProcStats(surfMesh, meshBb);
286     Info<< "Writing surface." << nl << endl;
287     surfMesh.searchableSurface::write();
289     Info<< "End\n" << endl;
291     return 0;
295 // ************************************************************************* //