histogram.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2016-2021 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 "histogram.H"
27 #include "coordSet.H"
28 #include "volFields.H"
30 
31 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35 namespace functionObjects
36 {
37  defineTypeNameAndDebug(histogram, 0);
38  addToRunTimeSelectionTable(functionObject, histogram, dictionary);
39 }
40 }
41 
42 
43 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
44 
46 (
47  const word& name,
48  const Time& runTime,
49  const dictionary& dict
50 )
51 :
52  fvMeshFunctionObject(name, runTime, dict),
53  file_(obr_, name)
54 {
55  read(dict);
56 }
57 
58 
59 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
60 
62 {}
63 
64 
65 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
66 
68 {
69  dict.lookup("field") >> fieldName_;
70  dict.lookup("max") >> max_;
71  min_ = dict.lookupOrDefault<scalar>("min", 0);
72  dict.lookup("nBins") >> nBins_;
73 
74  formatterPtr_ = setWriter::New(dict.lookup("setFormat"), dict);
75 
76  return true;
77 }
78 
79 
81 {
82  return wordList(fieldName_);
83 }
84 
85 
87 {
88  return true;
89 }
90 
91 
93 {
94  Log << type() << " " << name() << " write:" << nl;
95 
96  autoPtr<volScalarField> fieldPtr;
97  if (obr_.foundObject<volScalarField>(fieldName_))
98  {
99  Log << " Looking up field " << fieldName_ << endl;
100  }
101  else
102  {
103  Log << " Reading field " << fieldName_ << endl;
104  fieldPtr.reset
105  (
106  new volScalarField
107  (
108  IOobject
109  (
110  fieldName_,
111  mesh_.time().timeName(),
112  mesh_,
115  ),
116  mesh_
117  )
118  );
119  }
120 
121  const volScalarField& field =
122  (
123  fieldPtr.valid()
124  ? fieldPtr()
125  : obr_.lookupObject<volScalarField>(fieldName_)
126  );
127 
128  // Calculate the mid-points of bins for the graph axis
129  scalarField xBin(nBins_);
130  const scalar delta = (max_- min_)/nBins_;
131  scalar x = min_ + 0.5*delta;
132  forAll(xBin, i)
133  {
134  xBin[i] = x;
135  x += delta;
136  }
137 
138  scalarField volFrac(nBins_, 0);
139  const scalarField& V = mesh_.V();
140 
141  forAll(field, celli)
142  {
143  const label bini = (field[celli] - min_)/delta;
144  if (bini >= 0 && bini < nBins_)
145  {
146  volFrac[bini] += V[celli];
147  }
148  }
149 
151 
152  if (Pstream::master())
153  {
154  const scalar sumVol = sum(volFrac);
155 
156  if (sumVol > small)
157  {
158  volFrac /= sumVol;
159 
160  formatterPtr_().write
161  (
162  file_.baseTimeDir(),
163  typeName,
164  coordSet(true, fieldName_, xBin),
165  field.name(),
166  volFrac
167  );
168  }
169  }
170 
171  return true;
172 }
173 
174 
175 // ************************************************************************* //
static autoPtr< setWriter > New(const word &writeType, const IOstream::streamFormat writeFormat=IOstream::ASCII, const IOstream::compressionType writeCompression=IOstream::UNCOMPRESSED)
Select given write options.
Definition: setWriter.C:279
scalar delta
dictionary dict
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
virtual bool write()
Calculate the histogram and write.
Definition: histogram.C:92
const word & name() const
Return name.
Definition: IOobject.H:315
void reset(T *=nullptr)
If object pointer already set, delete object and set to given.
Definition: autoPtrI.H:114
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
addToRunTimeSelectionTable(functionObject, Qdot, dictionary)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
static bool master(const label communicator=0)
Am I the master process.
Definition: UPstream.H:423
static void listCombineGather(const List< commsStruct > &comms, List< T > &Value, const CombineOp &cop, const int tag, const label comm)
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:69
virtual bool execute()
Execute, currently does nothing.
Definition: histogram.C:86
Macros for easy insertion into run-time selection tables.
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
virtual wordList fields() const
Return the list of fields required.
Definition: histogram.C:80
bool read(const char *, int32_t &)
Definition: int32IO.C:85
virtual bool read(const dictionary &)
Read the histogram data.
Definition: histogram.C:67
histogram(const word &name, const Time &runTime, const dictionary &dict)
Construct from Time and dictionary.
Definition: histogram.C:46
Holds list of sampling positions.
Definition: coordSet.H:50
A class for handling words, derived from string.
Definition: word.H:59
bool valid() const
Return true if the autoPtr valid (ie, the pointer is set)
Definition: autoPtrI.H:83
static const char nl
Definition: Ostream.H:260
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
defineTypeNameAndDebug(Qdot, 0)
List< word > wordList
A List of words.
Definition: fileName.H:54
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
#define Log
Report write to Foam::Info if the local log switch is true.
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:52
Specialisation of Foam::functionObject for an Foam::fvMesh, providing a reference to the Foam::fvMesh...
rDeltaTY field()
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:98
Namespace for OpenFOAM.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:864