src/sampling/probes/probes.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 \*---------------------------------------------------------------------------*/
  26
  27 #include "probes.H"
  28 #include "volFields.H"
  29 #include "dictionary.H"
  30 #include "Time.H"
  31 #include "IOmanip.H"
  32
  33 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
  34
  35 namespace Foam
  36 {
  37     defineTypeNameAndDebug(probes, 0);
  38 }
  39
  40 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
  41
  42 void Foam::probes::findCells(const fvMesh& mesh)
  43 {
  44     if (cellList_.empty())
  45     {
  46         cellList_.setSize(probeLocations_.size());
  47
  48         forAll(probeLocations_, probeI)
  49         {
  50             cellList_[probeI] = mesh.findCell(probeLocations_[probeI]);
  51
  52             if (debug && cellList_[probeI] != -1)
  53             {
  54                 Pout<< "probes : found point " << probeLocations_[probeI]
  55                     << " in cell " << cellList_[probeI] << endl;
  56             }
  57         }
  58
  59
  60         // Check if all probes have been found.
  61         forAll(cellList_, probeI)
  62         {
  63             label cellI = cellList_[probeI];
  64
  65             // Check at least one processor with cell.
  66             reduce(cellI, maxOp<label>());
  67
  68             if (cellI == -1)
  69             {
  70                 if (Pstream::master())
  71                 {
  72                     WarningIn("probes::read()")
  73                         << "Did not find location " << probeLocations_[probeI]
  74                         << " in any cell. Skipping location." << endl;
  75                 }
  76             }
  77             else
  78             {
  79                 // Make sure location not on two domains.
  80                 if (cellList_[probeI] != -1 && cellList_[probeI] != cellI)
  81                 {
  82                     WarningIn("probes::read()")
  83                         << "Location " << probeLocations_[probeI]
  84                         << " seems to be on multiple domains:"
  85                         << " cell " << cellList_[probeI]
  86                         << " on my domain " << Pstream::myProcNo()
  87                         << " and cell " << cellI << " on some other domain."
  88                         << endl
  89                         << "This might happen if the probe location is on"
  90                         << " a processor patch. Change the location slightly"
  91                         << " to prevent this." << endl;
  92                 }
  93             }
  94         }
  95     }
  96 }
  97
  98
  99 bool Foam::probes::checkFieldTypes()
 100 {
 101     wordList fieldTypes(fieldNames_.size());
 102
 103     // check files for a particular time
 104     if (loadFromFiles_)
 105     {
 106         forAll(fieldNames_, fieldI)
 107         {
 108             IOobject io
 109             (
 110                 fieldNames_[fieldI],
 111                 obr_.time().timeName(),
 112                 refCast<const polyMesh>(obr_),
 113                 IOobject::MUST_READ,
 114                 IOobject::NO_WRITE,
 115                 false
 116             );
 117
 118             if (io.headerOk())
 119             {
 120                 fieldTypes[fieldI] = io.headerClassName();
 121             }
 122             else
 123             {
 124                 fieldTypes[fieldI] = "(notFound)";
 125             }
 126         }
 127     }
 128     else
 129     {
 130         // check objectRegistry
 131         forAll(fieldNames_, fieldI)
 132         {
 133             objectRegistry::const_iterator iter =
 134                 obr_.find(fieldNames_[fieldI]);
 135
 136             if (iter != obr_.end())
 137             {
 138                 fieldTypes[fieldI] = iter()->type();
 139             }
 140             else
 141             {
 142                 fieldTypes[fieldI] = "(notFound)";
 143             }
 144         }
 145     }
 146
 147
 148     label nFields = 0;
 149
 150     // classify fieldTypes
 151     nFields += countFields(scalarFields_, fieldTypes);
 152     nFields += countFields(vectorFields_, fieldTypes);
 153     nFields += countFields(sphericalTensorFields_, fieldTypes);
 154     nFields += countFields(symmTensorFields_, fieldTypes);
 155     nFields += countFields(tensorFields_, fieldTypes);
 156
 157     // concatenate all the lists into foundFields
 158     wordList foundFields(nFields);
 159
 160     label fieldI = 0;
 161     forAll(scalarFields_, i)
 162     {
 163         foundFields[fieldI++] = scalarFields_[i];
 164     }
 165     forAll(vectorFields_, i)
 166     {
 167         foundFields[fieldI++] = vectorFields_[i];
 168     }
 169     forAll(sphericalTensorFields_, i)
 170     {
 171         foundFields[fieldI++] = sphericalTensorFields_[i];
 172     }
 173     forAll(symmTensorFields_, i)
 174     {
 175         foundFields[fieldI++] = symmTensorFields_[i];
 176     }
 177     forAll(tensorFields_, i)
 178     {
 179         foundFields[fieldI++] = tensorFields_[i];
 180     }
 181
 182     if (Pstream::master())
 183     {
 184         fileName probeDir;
 185
 186         fileName probeSubDir = name_;
 187
 188         if (obr_.name() != polyMesh::defaultRegion)
 189         {
 190             probeSubDir = probeSubDir/obr_.name();
 191         }
 192         probeSubDir = probeSubDir/obr_.time().timeName();
 193
 194         if (Pstream::parRun())
 195         {
 196             // Put in undecomposed case
 197             // (Note: gives problems for distributed data running)
 198             probeDir = obr_.time().path()/".."/probeSubDir;
 199         }
 200         else
 201         {
 202             probeDir = obr_.time().path()/probeSubDir;
 203         }
 204
 205         // Close the file if any fields have been removed.
 206         forAllIter(HashPtrTable<OFstream>, probeFilePtrs_, iter)
 207         {
 208             if (findIndex(foundFields, iter.key()) == -1)
 209             {
 210                 if (debug)
 211                 {
 212                     Pout<< "close stream: " << iter()->name() << endl;
 213                 }
 214
 215                 delete probeFilePtrs_.remove(iter);
 216             }
 217         }
 218
 219         // Open new files for new fields. Keep existing files.
 220
 221         probeFilePtrs_.resize(2*foundFields.size());
 222
 223         forAll(foundFields, fieldI)
 224         {
 225             const word& fldName = foundFields[fieldI];
 226
 227             // Check if added field. If so open a stream for it.
 228
 229             if (!probeFilePtrs_.found(fldName))
 230             {
 231                 // Create directory if does not exist.
 232                 mkDir(probeDir);
 233
 234                 OFstream* sPtr = new OFstream(probeDir/fldName);
 235
 236                 if (debug)
 237                 {
 238                     Pout<< "open  stream: " << sPtr->name() << endl;
 239                 }
 240
 241                 probeFilePtrs_.insert(fldName, sPtr);
 242
 243                 unsigned int w = IOstream::defaultPrecision() + 7;
 244
 245                 for (direction cmpt=0; cmpt<vector::nComponents; cmpt++)
 246                 {
 247                     *sPtr<< '#' << setw(IOstream::defaultPrecision() + 6)
 248                         << vector::componentNames[cmpt];
 249
 250                     forAll(probeLocations_, probeI)
 251                     {
 252                         *sPtr<< ' ' << setw(w) << probeLocations_[probeI][cmpt];
 253                     }
 254                     *sPtr << endl;
 255                 }
 256
 257                 *sPtr<< '#' << setw(IOstream::defaultPrecision() + 6)
 258                     << "Time" << endl;
 259             }
 260         }
 261
 262         if (debug)
 263         {
 264             Pout<< "Probing fields:" << foundFields << nl
 265                 << "Probing locations:" << probeLocations_ << nl
 266                 << endl;
 267         }
 268     }
 269
 270
 271     return nFields > 0;
 272 }
 273
 274
 275 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 276
 277 Foam::probes::probes
 278 (
 279     const word& name,
 280     const objectRegistry& obr,
 281     const dictionary& dict,
 282     const bool loadFromFiles
 283 )
 284 :
 285     name_(name),
 286     obr_(obr),
 287     loadFromFiles_(loadFromFiles),
 288     fieldNames_(),
 289     probeLocations_(),
 290     scalarFields_(),
 291     vectorFields_(),
 292     sphericalTensorFields_(),
 293     symmTensorFields_(),
 294     tensorFields_(),
 295     cellList_(),
 296     probeFilePtrs_()
 297 {
 298     read(dict);
 299 }
 300
 301
 302 // * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
 303
 304 Foam::probes::~probes()
 305 {}
 306
 307
 308 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 309
 310 void Foam::probes::execute()
 311 {
 312     // Do nothing - only valid on write
 313 }
 314
 315
 316 void Foam::probes::end()
 317 {
 318     // Do nothing - only valid on write
 319 }
 320
 321
 322 void Foam::probes::write()
 323 {
 324     if (probeLocations_.size() && checkFieldTypes())
 325     {
 326         sampleAndWrite(scalarFields_);
 327         sampleAndWrite(vectorFields_);
 328         sampleAndWrite(sphericalTensorFields_);
 329         sampleAndWrite(symmTensorFields_);
 330         sampleAndWrite(tensorFields_);
 331     }
 332 }
 333
 334
 335 void Foam::probes::read(const dictionary& dict)
 336 {
 337     dict.lookup("fields") >> fieldNames_;
 338     dict.lookup("probeLocations") >> probeLocations_;
 339
 340     // Force all cell locations to be redetermined
 341     cellList_.clear();
 342     findCells(refCast<const fvMesh>(obr_));
 343     checkFieldTypes();
 344 }
 345
 346
 347 // ************************************************************************* //