src/sampling/probes/probes.C

   1 /*---------------------------------------------------------------------------*\
   2   =========                 |
   3   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   4    \\    /   O peration     |
   5     \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
   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 "probes.H"
  27 #include "volFields.H"
  28 #include "dictionary.H"
  29 #include "Time.H"
  30 #include "IOmanip.H"
  31
  32 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
  33
  34 defineTypeNameAndDebug(Foam::probes, 0);
  35
  36
  37 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
  38
  39 void Foam::probes::findElements(const fvMesh& mesh)
  40 {
  41     elementList_.clear();
  42     elementList_.setSize(size());
  43
  44     forAll(*this, probeI)
  45     {
  46         const vector& location = operator[](probeI);
  47
  48         elementList_[probeI] = mesh.findCell(location);
  49
  50         if (debug && elementList_[probeI] != -1)
  51         {
  52             Pout<< "probes : found point " << location
  53                 << " in cell " << elementList_[probeI] << endl;
  54         }
  55     }
  56
  57
  58     // Check if all probes have been found.
  59     forAll(elementList_, probeI)
  60     {
  61         const vector& location = operator[](probeI);
  62         label cellI = elementList_[probeI];
  63
  64         // Check at least one processor with cell.
  65         reduce(cellI, maxOp<label>());
  66
  67         if (cellI == -1)
  68         {
  69             if (Pstream::master())
  70             {
  71                 WarningIn("probes::read()")
  72                     << "Did not find location " << location
  73                     << " in any cell. Skipping location." << endl;
  74             }
  75         }
  76         else
  77         {
  78             // Make sure location not on two domains.
  79             if (elementList_[probeI] != -1 && elementList_[probeI] != cellI)
  80             {
  81                 WarningIn("probes::read()")
  82                     << "Location " << location
  83                     << " seems to be on multiple domains:"
  84                     << " cell " << elementList_[probeI]
  85                     << " on my domain " << Pstream::myProcNo()
  86                         << " and cell " << cellI << " on some other domain."
  87                     << endl
  88                     << "This might happen if the probe location is on"
  89                     << " a processor patch. Change the location slightly"
  90                     << " to prevent this." << endl;
  91             }
  92         }
  93     }
  94 }
  95
  96
  97 Foam::label Foam::probes::prepare()
  98 {
  99     const label nFields = classifyFields();
 100
 101     // adjust file streams
 102     if (Pstream::master())
 103     {
 104         wordHashSet currentFields;
 105
 106         currentFields.insert(scalarFields_);
 107         currentFields.insert(vectorFields_);
 108         currentFields.insert(sphericalTensorFields_);
 109         currentFields.insert(symmTensorFields_);
 110         currentFields.insert(tensorFields_);
 111
 112         if (debug)
 113         {
 114             Info<< "Probing fields:" << currentFields << nl
 115                 << "Probing locations:" << *this << nl
 116                 << endl;
 117         }
 118
 119
 120         fileName probeDir;
 121         fileName probeSubDir = name_;
 122
 123         if (mesh_.name() != polyMesh::defaultRegion)
 124         {
 125             probeSubDir = probeSubDir/mesh_.name();
 126         }
 127         probeSubDir = probeSubDir/mesh_.time().timeName();
 128
 129         if (Pstream::parRun())
 130         {
 131             // Put in undecomposed case
 132             // (Note: gives problems for distributed data running)
 133             probeDir = mesh_.time().path()/".."/probeSubDir;
 134         }
 135         else
 136         {
 137             probeDir = mesh_.time().path()/probeSubDir;
 138         }
 139
 140         // ignore known fields, close streams for fields that no longer exist
 141         forAllIter(HashPtrTable<OFstream>, probeFilePtrs_, iter)
 142         {
 143             if (!currentFields.erase(iter.key()))
 144             {
 145                 if (debug)
 146                 {
 147                     Info<< "close probe stream: " << iter()->name() << endl;
 148                 }
 149
 150                 delete probeFilePtrs_.remove(iter);
 151             }
 152         }
 153
 154         // currentFields now just has the new fields - open streams for them
 155         forAllConstIter(wordHashSet, currentFields, iter)
 156         {
 157             const word& fieldName = iter.key();
 158
 159             // Create directory if does not exist.
 160             mkDir(probeDir);
 161
 162             OFstream* sPtr = new OFstream(probeDir/fieldName);
 163
 164             if (debug)
 165             {
 166                 Info<< "open probe stream: " << sPtr->name() << endl;
 167             }
 168
 169             probeFilePtrs_.insert(fieldName, sPtr);
 170
 171             unsigned int w = IOstream::defaultPrecision() + 7;
 172
 173             for (direction cmpt=0; cmpt<vector::nComponents; cmpt++)
 174             {
 175                 *sPtr<< '#' << setw(IOstream::defaultPrecision() + 6)
 176                     << vector::componentNames[cmpt];
 177
 178                 forAll(*this, probeI)
 179                 {
 180                     *sPtr<< ' ' << setw(w) << operator[](probeI)[cmpt];
 181                 }
 182                 *sPtr << endl;
 183             }
 184
 185             *sPtr<< '#' << setw(IOstream::defaultPrecision() + 6)
 186                 << "Time" << endl;
 187         }
 188     }
 189
 190     return nFields;
 191 }
 192
 193
 194 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 195
 196 Foam::probes::probes
 197 (
 198     const word& name,
 199     const objectRegistry& obr,
 200     const dictionary& dict,
 201     const bool loadFromFiles
 202 )
 203 :
 204     pointField(0),
 205     name_(name),
 206     mesh_(refCast<const fvMesh>(obr)),
 207     loadFromFiles_(loadFromFiles)
 208 {
 209     read(dict);
 210 }
 211
 212
 213 // * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
 214
 215 Foam::probes::~probes()
 216 {}
 217
 218
 219 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 220
 221 void Foam::probes::execute()
 222 {
 223     // Do nothing - only valid on write
 224 }
 225
 226
 227 void Foam::probes::end()
 228 {
 229     // Do nothing - only valid on write
 230 }
 231
 232
 233 void Foam::probes::write()
 234 {
 235     if (size() && prepare())
 236     {
 237         sampleAndWrite(scalarFields_);
 238         sampleAndWrite(vectorFields_);
 239         sampleAndWrite(sphericalTensorFields_);
 240         sampleAndWrite(symmTensorFields_);
 241         sampleAndWrite(tensorFields_);
 242     }
 243 }
 244
 245
 246 void Foam::probes::read(const dictionary& dict)
 247 {
 248     dict.lookup("probeLocations") >> *this;
 249     dict.lookup("fields") >> fieldSelection_;
 250
 251     // redetermined all cell locations
 252     findElements(mesh_);
 253     prepare();
 254 }
 255
 256
 257 // ************************************************************************* //