src/sampling/probes/probes.C

   1 /*---------------------------------------------------------------------------*\
   2   =========                 |
   3   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   4    \\    /   O peration     |
   5     \\  /    A nd           | Copyright (C) 1991-2008 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 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_.size() == 0)
  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         if (Pstream::parRun())
 186         {
 187             // Put in undecomposed case
 188             // (Note: gives problems for distributed data running)
 189             probeDir = obr_.time().path()/".."/name_/obr_.time().timeName();
 190         }
 191         else
 192         {
 193             probeDir = obr_.time().path()/name_/obr_.time().timeName();
 194         }
 195
 196         // Close the file if any fields have been removed.
 197         forAllIter(HashPtrTable<OFstream>, probeFilePtrs_, iter)
 198         {
 199             if (findIndex(foundFields, iter.key()) == -1)
 200             {
 201                 if (debug)
 202                 {
 203                     Pout<< "close stream: " << iter()->name() << endl;
 204                 }
 205
 206                 delete probeFilePtrs_.remove(iter);
 207             }
 208         }
 209
 210         // Open new files for new fields. Keep existing files.
 211
 212         probeFilePtrs_.resize(2*foundFields.size());
 213
 214         forAll(foundFields, fieldI)
 215         {
 216             const word& fldName = foundFields[fieldI];
 217
 218             // Check if added field. If so open a stream for it.
 219
 220             if (!probeFilePtrs_.found(fldName))
 221             {
 222                 // Create directory if does not exist.
 223                 mkDir(probeDir);
 224
 225                 OFstream* sPtr = new OFstream(probeDir/fldName);
 226
 227                 if (debug)
 228                 {
 229                     Pout<< "open  stream: " << sPtr->name() << endl;
 230                 }
 231
 232                 probeFilePtrs_.insert(fldName, sPtr);
 233
 234                 unsigned int w = IOstream::defaultPrecision() + 7;
 235
 236                 for (direction cmpt=0; cmpt<vector::nComponents; cmpt++)
 237                 {
 238                     *sPtr<< '#' << setw(IOstream::defaultPrecision() + 6)
 239                         << vector::componentNames[cmpt];
 240
 241                     forAll(probeLocations_, probeI)
 242                     {
 243                         *sPtr<< ' ' << setw(w) << probeLocations_[probeI][cmpt];
 244                     }
 245                     *sPtr << endl;
 246                 }
 247
 248                 *sPtr<< '#' << setw(IOstream::defaultPrecision() + 6)
 249                     << "Time" << endl;
 250             }
 251         }
 252
 253         if (debug)
 254         {
 255             Pout<< "Probing fields:" << foundFields << nl
 256                 << "Probing locations:" << probeLocations_ << nl
 257                 << endl;
 258         }
 259     }
 260
 261
 262     return nFields > 0;
 263 }
 264
 265
 266 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 267
 268 Foam::probes::probes
 269 (
 270     const word& name,
 271     const objectRegistry& obr,
 272     const dictionary& dict,
 273     const bool loadFromFiles
 274 )
 275 :
 276     name_(name),
 277     obr_(obr),
 278     loadFromFiles_(loadFromFiles),
 279     fieldNames_(0),
 280     probeLocations_(0),
 281     scalarFields_(),
 282     vectorFields_(),
 283     sphericalTensorFields_(),
 284     symmTensorFields_(),
 285     tensorFields_(),
 286     cellList_(0),
 287     probeFilePtrs_(0)
 288 {
 289     read(dict);
 290 }
 291
 292
 293 // * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
 294
 295 Foam::probes::~probes()
 296 {}
 297
 298
 299 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 300
 301 void Foam::probes::execute()
 302 {
 303     // Do nothing - only valid on write
 304 }
 305
 306
 307 void Foam::probes::write()
 308 {
 309     if (probeLocations_.size() && checkFieldTypes())
 310     {
 311         sampleAndWrite(scalarFields_);
 312         sampleAndWrite(vectorFields_);
 313         sampleAndWrite(sphericalTensorFields_);
 314         sampleAndWrite(symmTensorFields_);
 315         sampleAndWrite(tensorFields_);
 316     }
 317 }
 318
 319
 320 void Foam::probes::read(const dictionary& dict)
 321 {
 322     dict.lookup("fields") >> fieldNames_;
 323     dict.lookup("probeLocations") >> probeLocations_;
 324
 325     // Force all cell locations to be redetermined
 326     cellList_.clear();
 327     findCells(refCast<const fvMesh>(obr_));
 328     checkFieldTypes();
 329 }
 330
 331
 332 // ************************************************************************* //