src/sampling/sampledSet/sampledSets/sampledSets.C

   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 \*---------------------------------------------------------------------------*/
  25
  26 #include "sampledSets.H"
  27 #include "dictionary.H"
  28 #include "Time.H"
  29 #include "volFields.H"
  30 #include "ListListOps.H"
  31 #include "SortableList.H"
  32 #include "volPointInterpolation.H"
  33
  34 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
  35
  36 namespace Foam
  37 {
  38     defineTypeNameAndDebug(sampledSets, 0);
  39 }
  40
  41 bool Foam::sampledSets::verbose_ = false;
  42
  43
  44 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
  45
  46 bool Foam::sampledSets::checkFieldTypes()
  47 {
  48     wordList fieldTypes(fieldNames_.size());
  49
  50     // check files for a particular time
  51     if (loadFromFiles_)
  52     {
  53         forAll(fieldNames_, fieldi)
  54         {
  55             IOobject io
  56             (
  57                 fieldNames_[fieldi],
  58                 mesh_.time().timeName(),
  59                 mesh_,
  60                 IOobject::MUST_READ,
  61                 IOobject::NO_WRITE,
  62                 false
  63             );
  64
  65             if (io.headerOk())
  66             {
  67                 fieldTypes[fieldi] = io.headerClassName();
  68             }
  69             else
  70             {
  71                 fieldTypes[fieldi] = "(notFound)";
  72             }
  73         }
  74     }
  75     else
  76     {
  77         // check objectRegistry
  78         forAll(fieldNames_, fieldi)
  79         {
  80             objectRegistry::const_iterator iter =
  81                 mesh_.find(fieldNames_[fieldi]);
  82
  83             if (iter != mesh_.objectRegistry::end())
  84             {
  85                 fieldTypes[fieldi] = iter()->type();
  86             }
  87             else
  88             {
  89                 fieldTypes[fieldi] = "(notFound)";
  90             }
  91         }
  92     }
  93
  94
  95     label nFields = 0;
  96
  97     // classify fieldTypes
  98     nFields += grep(scalarFields_, fieldTypes);
  99     nFields += grep(vectorFields_, fieldTypes);
 100     nFields += grep(sphericalTensorFields_, fieldTypes);
 101     nFields += grep(symmTensorFields_, fieldTypes);
 102     nFields += grep(tensorFields_, fieldTypes);
 103
 104     if (Pstream::master())
 105     {
 106         if (debug)
 107         {
 108             Pout<< "timeName = " << mesh_.time().timeName() << nl
 109                 << "scalarFields    " << scalarFields_ << nl
 110                 << "vectorFields    " << vectorFields_ << nl
 111                 << "sphTensorFields " << sphericalTensorFields_ << nl
 112                 << "symTensorFields " << symmTensorFields_ <<nl
 113                 << "tensorFields    " << tensorFields_ <<nl;
 114         }
 115
 116         if (nFields > 0)
 117         {
 118             if (debug)
 119             {
 120                 Pout<< "Creating directory "
 121                     << outputPath_/mesh_.time().timeName()
 122                     << nl << endl;
 123             }
 124
 125             mkDir(outputPath_/mesh_.time().timeName());
 126         }
 127     }
 128
 129     return nFields > 0;
 130 }
 131
 132
 133 void Foam::sampledSets::combineSampledSets
 134 (
 135     PtrList<coordSet>& masterSampledSets,
 136     labelListList& indexSets
 137 )
 138 {
 139     // Combine sampleSets from processors. Sort by curveDist. Return
 140     // ordering in indexSets.
 141     // Note: only master results are valid
 142
 143     masterSampledSets_.clear();
 144     masterSampledSets_.setSize(size());
 145     indexSets_.setSize(size());
 146
 147     const PtrList<sampledSet>& sampledSets = *this;
 148
 149     forAll(sampledSets, seti)
 150     {
 151         const sampledSet& samplePts = sampledSets[seti];
 152
 153         // Collect data from all processors
 154         List<List<point> > gatheredPts(Pstream::nProcs());
 155         gatheredPts[Pstream::myProcNo()] = samplePts;
 156         Pstream::gatherList(gatheredPts);
 157
 158         List<labelList> gatheredSegments(Pstream::nProcs());
 159         gatheredSegments[Pstream::myProcNo()] = samplePts.segments();
 160         Pstream::gatherList(gatheredSegments);
 161
 162         List<scalarList> gatheredDist(Pstream::nProcs());
 163         gatheredDist[Pstream::myProcNo()] = samplePts.curveDist();
 164         Pstream::gatherList(gatheredDist);
 165
 166
 167         // Combine processor lists into one big list.
 168         List<point> allPts
 169         (
 170             ListListOps::combine<List<point> >
 171             (
 172                 gatheredPts, accessOp<List<point> >()
 173             )
 174         );
 175         labelList allSegments
 176         (
 177             ListListOps::combine<labelList>
 178             (
 179                 gatheredSegments, accessOp<labelList>()
 180             )
 181         );
 182         scalarList allCurveDist
 183         (
 184             ListListOps::combine<scalarList>
 185             (
 186                 gatheredDist, accessOp<scalarList>()
 187             )
 188         );
 189
 190
 191         if (Pstream::master() && allCurveDist.size() == 0)
 192         {
 193             WarningIn("sampledSets::combineSampledSets(..)")
 194                 << "Sample set " << samplePts.name()
 195                 << " has zero points." << endl;
 196         }
 197
 198         // Sort curveDist and use to fill masterSamplePts
 199         SortableList<scalar> sortedDist(allCurveDist);
 200         indexSets[seti] = sortedDist.indices();
 201
 202         // Get reference point (note: only master has all points)
 203         point refPt;
 204
 205         if (allPts.size())
 206         {
 207             refPt = samplePts.getRefPoint(allPts);
 208         }
 209         else
 210         {
 211             refPt = vector::zero;
 212         }
 213
 214
 215         masterSampledSets.set
 216         (
 217             seti,
 218             new coordSet
 219             (
 220                 samplePts.name(),
 221                 samplePts.axis(),
 222                 UIndirectList<point>(allPts, indexSets[seti]),
 223                 refPt
 224             )
 225         );
 226     }
 227 }
 228
 229
 230 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 231
 232 Foam::sampledSets::sampledSets
 233 (
 234     const word& name,
 235     const objectRegistry& obr,
 236     const dictionary& dict,
 237     const bool loadFromFiles
 238 )
 239 :
 240     PtrList<sampledSet>(),
 241     name_(name),
 242     mesh_(refCast<const fvMesh>(obr)),
 243     loadFromFiles_(loadFromFiles),
 244     outputPath_(fileName::null),
 245     searchEngine_(mesh_, true),
 246     fieldNames_(),
 247     interpolationScheme_(word::null),
 248     writeFormat_(word::null)
 249 {
 250     if (Pstream::parRun())
 251     {
 252         outputPath_ = mesh_.time().path()/".."/name_;
 253     }
 254     else
 255     {
 256         outputPath_ = mesh_.time().path()/name_;
 257     }
 258     if (mesh_.name() != fvMesh::defaultRegion)
 259     {
 260         outputPath_ = outputPath_/mesh_.name();
 261     }
 262
 263     read(dict);
 264 }
 265
 266
 267 // * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
 268
 269 Foam::sampledSets::~sampledSets()
 270 {}
 271
 272
 273 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 274
 275 void Foam::sampledSets::verbose(const bool verbosity)
 276 {
 277     verbose_ = verbosity;
 278 }
 279
 280
 281 void Foam::sampledSets::execute()
 282 {
 283     // Do nothing - only valid on write
 284 }
 285
 286
 287 void Foam::sampledSets::end()
 288 {
 289     // Do nothing - only valid on write
 290 }
 291
 292
 293 void Foam::sampledSets::write()
 294 {
 295     if (size() && checkFieldTypes())
 296     {
 297         sampleAndWrite(scalarFields_);
 298         sampleAndWrite(vectorFields_);
 299         sampleAndWrite(sphericalTensorFields_);
 300         sampleAndWrite(symmTensorFields_);
 301         sampleAndWrite(tensorFields_);
 302     }
 303 }
 304
 305
 306 void Foam::sampledSets::read(const dictionary& dict)
 307 {
 308     dict_ = dict;
 309
 310     fieldNames_ = wordList(dict_.lookup("fields"));
 311
 312     interpolationScheme_ = "cell";
 313     dict_.readIfPresent("interpolationScheme", interpolationScheme_);
 314
 315     dict_.lookup("setFormat") >> writeFormat_;
 316
 317     scalarFields_.clear();
 318     vectorFields_.clear();
 319     sphericalTensorFields_.clear();
 320     symmTensorFields_.clear();
 321     tensorFields_.clear();
 322
 323     PtrList<sampledSet> newList
 324     (
 325         dict_.lookup("sets"),
 326         sampledSet::iNew(mesh_, searchEngine_)
 327     );
 328     transfer(newList);
 329     combineSampledSets(masterSampledSets_, indexSets_);
 330
 331     if (Pstream::master() && debug)
 332     {
 333         Pout<< "sample fields:" << fieldNames_ << nl
 334             << "sample sets:" << nl << "(" << nl;
 335
 336         forAll(*this, si)
 337         {
 338             Pout << "  " << operator[](si) << endl;
 339         }
 340         Pout << ")" << endl;
 341     }
 342 }
 343
 344
 345 void Foam::sampledSets::correct()
 346 {
 347     // reset interpolation
 348     pointMesh::Delete(mesh_);
 349     volPointInterpolation::Delete(mesh_);
 350
 351     searchEngine_.correct();
 352
 353     PtrList<sampledSet> newList
 354     (
 355         dict_.lookup("sets"),
 356         sampledSet::iNew(mesh_, searchEngine_)
 357     );
 358     transfer(newList);
 359     combineSampledSets(masterSampledSets_, indexSets_);
 360 }
 361
 362
 363 void Foam::sampledSets::updateMesh(const mapPolyMesh&)
 364 {
 365     correct();
 366 }
 367
 368
 369 void Foam::sampledSets::movePoints(const pointField&)
 370 {
 371     correct();
 372 }
 373
 374
 375 void Foam::sampledSets::readUpdate(const polyMesh::readUpdateState state)
 376 {
 377     if (state != polyMesh::UNCHANGED)
 378     {
 379         correct();
 380     }
 381 }
 382
 383
 384 // ************************************************************************* //