applications/utilities/surface/surfaceRedistributePar/surfaceRedistributePar.C

   1 /*---------------------------------------------------------------------------*\
   2   =========                 |
   3   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   4    \\    /   O peration     |
   5     \\  /    A nd           | Copyright held by original author
   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 the
  13     Free Software Foundation; either version 2 of the License, or (at your
  14     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, write to the Free Software Foundation,
  23     Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
  24
  25 Application
  26     surfaceRedistributePar
  27
  28 Description
  29     (Re)distribution of triSurface. Either takes an undecomposed surface
  30     or an already decomposed surface and redistribute it so each processor
  31     has all triangles that overlap its mesh.
  32
  33 Note
  34     - best decomposition option is hierarchGeomDecomp since
  35     guarantees square decompositions.
  36     - triangles might be present on multiple processors.
  37     - merging uses geometric tolerance so take care with writing precision.
  38
  39 \*---------------------------------------------------------------------------*/
  40
  41 #include "treeBoundBox.H"
  42 #include "FixedList.H"
  43 #include "argList.H"
  44 #include "objectRegistry.H"
  45 #include "Time.H"
  46 #include "polyMesh.H"
  47 #include "distributedTriSurfaceMesh.H"
  48 #include "mapDistribute.H"
  49 #include "triSurfaceFields.H"
  50 #include "Pair.H"
  51
  52 using namespace Foam;
  53
  54 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
  55
  56 // Print on master all the per-processor surface stats.
  57 void writeProcStats
  58 (
  59     const triSurface& s,
  60     const List<List<treeBoundBox> >& meshBb
  61 )
  62 {
  63     // Determine surface bounding boxes, faces, points
  64     List<treeBoundBox> surfBb(Pstream::nProcs());
  65     {
  66         surfBb[Pstream::myProcNo()] = treeBoundBox(s.points());
  67         Pstream::gatherList(surfBb);
  68         Pstream::scatterList(surfBb);
  69     }
  70
  71     labelList nPoints(Pstream::nProcs());
  72     nPoints[Pstream::myProcNo()] = s.points().size();
  73     Pstream::gatherList(nPoints);
  74     Pstream::scatterList(nPoints);
  75
  76     labelList nFaces(Pstream::nProcs());
  77     nFaces[Pstream::myProcNo()] = s.size();
  78     Pstream::gatherList(nFaces);
  79     Pstream::scatterList(nFaces);
  80
  81     forAll(surfBb, procI)
  82     {
  83         const List<treeBoundBox>& bbs = meshBb[procI];
  84
  85         Info<< "processor" << procI << endl
  86             << "\tMesh bounds          : " << bbs[0] << nl;
  87         for (label i = 1; i < bbs.size(); i++)
  88         {
  89             Info<< "\t                       " << bbs[i]<< nl;
  90         }
  91         Info<< "\tSurface bounding box : " << surfBb[procI] << nl
  92             << "\tTriangles            : " << nFaces[procI] << nl
  93             << "\tVertices             : " << nPoints[procI]
  94             << endl;
  95     }
  96     Info<< endl;
  97 }
  98
  99
 100 // Main program:
 101
 102 int main(int argc, char *argv[])
 103 {
 104     argList::validArgs.append("triSurfaceMesh");
 105     argList::validArgs.append("distributionType");
 106
 107     argList::validOptions.insert("keepNonMapped", "");
 108 #   include "setRootCase.H"
 109 #   include "createTime.H"
 110     runTime.functionObjects().off();
 111
 112     fileName surfFileName(args.additionalArgs()[0]);
 113     Info<< "Reading surface from " << surfFileName << nl << endl;
 114
 115     const word distType(args.additionalArgs()[1]);
 116
 117     Info<< "Using distribution method "
 118         << distributedTriSurfaceMesh::distributionTypeNames_[distType]
 119         << " " << distType << nl << endl;
 120
 121     bool keepNonMapped = args.options().found("keepNonMapped");
 122
 123     if (keepNonMapped)
 124     {
 125         Info<< "Preserving surface outside of mesh bounds." << nl << endl;
 126     }
 127     else
 128     {
 129         Info<< "Removing surface outside of mesh bounds." << nl << endl;
 130     }
 131
 132
 133     if (!Pstream::parRun())
 134     {
 135         FatalErrorIn(args.executable())
 136             << "Please run this program on the decomposed case."
 137             << " It will read surface " << surfFileName
 138             << " and decompose it such that it overlaps the mesh bounding box."
 139             << exit(FatalError);
 140     }
 141
 142
 143 #   include "createPolyMesh.H"
 144
 145     Random rndGen(653213);
 146
 147     // Determine mesh bounding boxes:
 148     List<List<treeBoundBox> > meshBb(Pstream::nProcs());
 149     {
 150         meshBb[Pstream::myProcNo()] = List<treeBoundBox>
 151         (
 152             1,
 153             treeBoundBox
 154             (
 155                 boundBox(mesh.points(), false)
 156             ).extend(rndGen, 1E-3)
 157         );
 158         Pstream::gatherList(meshBb);
 159         Pstream::scatterList(meshBb);
 160     }
 161
 162     IOobject io
 163     (
 164         surfFileName,         // name
 165         //runTime.findInstance("triSurface", surfFileName),   // instance
 166         runTime.constant(),   // instance
 167         "triSurface",         // local
 168         runTime,              // registry
 169         IOobject::MUST_READ,
 170         IOobject::NO_WRITE
 171     );
 172
 173     const fileName actualPath(io.filePath());
 174     fileName localPath(actualPath);
 175     localPath.replace(runTime.rootPath() + '/', "");
 176
 177     if (actualPath == io.objectPath())
 178     {
 179         Info<< "Loading local (decomposed) surface " << localPath << nl <<endl;
 180     }
 181     else
 182     {
 183         Info<< "Loading undecomposed surface " << localPath << nl << endl;
 184     }
 185
 186
 187     // Create dummy dictionary for bounding boxes if does not exist.
 188     if (!isFile(actualPath / "Dict"))
 189     {
 190         dictionary dict;
 191         dict.add("bounds", meshBb[Pstream::myProcNo()]);
 192         dict.add("distributionType", distType);
 193         dict.add("mergeDistance", SMALL);
 194
 195         IOdictionary ioDict
 196         (
 197             IOobject
 198             (
 199                 io.name() + "Dict",
 200                 io.instance(),
 201                 io.local(),
 202                 io.db(),
 203                 IOobject::NO_READ,
 204                 IOobject::NO_WRITE,
 205                 false
 206             ),
 207             dict
 208         );
 209
 210         Info<< "Writing dummy bounds dictionary to " << ioDict.name()
 211             << nl << endl;
 212
 213         ioDict.regIOobject::writeObject
 214         (
 215             IOstream::ASCII,
 216             IOstream::currentVersion,
 217             ioDict.time().writeCompression()
 218         );
 219     }
 220
 221
 222     // Load surface
 223     distributedTriSurfaceMesh surfMesh(io);
 224     Info<< "Loaded surface" << nl << endl;
 225
 226
 227     // Generate a test field
 228     {
 229         const triSurface& s = static_cast<const triSurface&>(surfMesh);
 230
 231         autoPtr<triSurfaceVectorField> fcPtr
 232         (
 233             new triSurfaceVectorField
 234             (
 235                 IOobject
 236                 (
 237                     surfMesh.searchableSurface::name(),     // name
 238                     surfMesh.searchableSurface::instance(), // instance
 239                     surfMesh.searchableSurface::local(),    // local
 240                     surfMesh,
 241                     IOobject::NO_READ,
 242                     IOobject::AUTO_WRITE
 243                 ),
 244                 surfMesh,
 245                 dimLength
 246             )
 247         );
 248         triSurfaceVectorField& fc = fcPtr();
 249
 250         forAll(fc, triI)
 251         {
 252             fc[triI] = s[triI].centre(s.points());
 253         }
 254
 255         // Steal pointer and store object on surfMesh
 256         fcPtr.ptr()->store();
 257     }
 258
 259
 260     // Write per-processor stats
 261     Info<< "Before redistribution:" << endl;
 262     writeProcStats(surfMesh, meshBb);
 263
 264
 265     // Do redistribution
 266     Info<< "Redistributing surface" << nl << endl;
 267     autoPtr<mapDistribute> faceMap;
 268     autoPtr<mapDistribute> pointMap;
 269     surfMesh.distribute
 270     (
 271         meshBb[Pstream::myProcNo()],
 272         keepNonMapped,
 273         faceMap,
 274         pointMap
 275     );
 276     faceMap.clear();
 277     pointMap.clear();
 278
 279     Info<< endl;
 280
 281
 282     // Write per-processor stats
 283     Info<< "After redistribution:" << endl;
 284     writeProcStats(surfMesh, meshBb);
 285
 286
 287     Info<< "Writing surface." << nl << endl;
 288     surfMesh.searchableSurface::write();
 289
 290     Info<< "End\n" << endl;
 291
 292     return 0;
 293 }
 294
 295
 296 // ************************************************************************* //