histogram.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 2016 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 "volFields.H"
29 
30 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
31 
32 namespace Foam
33 {
34 namespace functionObjects
35 {
36  defineTypeNameAndDebug(histogram, 0);
37  addToRunTimeSelectionTable(functionObject, histogram, dictionary);
38 }
39 }
40 
41 
42 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
43 
44 void Foam::functionObjects::histogram::writeGraph
45 (
46  const coordSet& coords,
47  const word& fieldName,
48  const scalarField& values
49 ) const
50 {
51  const wordList fieldNames(1, fieldName);
52 
53  fileName outputPath = baseTimeDir();
54  mkDir(outputPath);
55  OFstream graphFile
56  (
57  outputPath/formatterPtr_().getFileName(coords, fieldNames)
58  );
59 
60  Log << " Writing histogram of " << fieldName
61  << " to " << graphFile.name() << endl;
62 
63  List<const scalarField*> yPtrs(1);
64  yPtrs[0] = &values;
65  formatterPtr_().write(coords, fieldNames, yPtrs, graphFile);
66 }
67 
68 
69 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
70 
71 Foam::functionObjects::histogram::histogram
72 (
73  const word& name,
74  const Time& runTime,
75  const dictionary& dict
76 )
77 :
78  writeFile(name, runTime, dict, name)
79 {
80  if (!isA<fvMesh>(obr_))
81  {
83  << "objectRegistry is not an fvMesh" << exit(FatalError);
84  }
85 
86  read(dict);
87 }
88 
89 
90 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
91 
93 {}
94 
95 
96 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
97 
99 {
100  dict.lookup("field") >> fieldName_;
101  dict.lookup("max") >> max_;
102  min_ = dict.lookupOrDefault<scalar>("min", 0);
103  dict.lookup("nBins") >> nBins_;
104 
105  word format(dict.lookup("setFormat"));
106  formatterPtr_ = writer<scalar>::New(format);
107 
108  return true;
109 }
110 
111 
113 {
114  return true;
115 }
116 
117 
119 {
120  Log << type() << " " << name() << " write:" << nl;
121 
122  const fvMesh& mesh = refCast<const fvMesh>(obr_);
123 
124  autoPtr<volScalarField> fieldPtr;
125  if (obr_.foundObject<volScalarField>(fieldName_))
126  {
127  Log << " Looking up field " << fieldName_ << endl;
128  }
129  else
130  {
131  Log << " Reading field " << fieldName_ << endl;
132  fieldPtr.reset
133  (
134  new volScalarField
135  (
136  IOobject
137  (
138  fieldName_,
139  mesh.time().timeName(),
140  mesh,
143  ),
144  mesh
145  )
146  );
147  }
148 
149  const volScalarField& field =
150  (
151  fieldPtr.valid()
152  ? fieldPtr()
153  : obr_.lookupObject<volScalarField>(fieldName_)
154  );
155 
156  // Calculate the mid-points of bins for the graph axis
157  pointField xBin(nBins_);
158  const scalar delta = (max_- min_)/nBins_;
159 
160  scalar x = min_ + 0.5*delta;
161  forAll(xBin, i)
162  {
163  xBin[i] = point(x, 0, 0);
164  x += delta;
165  }
166 
167  scalarField volFrac(nBins_, 0);
168  const scalarField& V = mesh.V();
169 
170  forAll(field, celli)
171  {
172  const label bini = (field[celli] - min_)/delta;
173  if (bini >= 0 && bini < nBins_)
174  {
175  volFrac[bini] += V[celli];
176  }
177  }
178 
180 
181  if (Pstream::master())
182  {
183  const scalar sumVol = sum(volFrac);
184 
185  if (sumVol > SMALL)
186  {
187  volFrac /= sumVol;
188 
189  const coordSet coords
190  (
191  "Volume_Fraction",
192  "x",
193  xBin,
194  mag(xBin)
195  );
196 
197  writeGraph(coords, field.name(), volFrac);
198  }
199  }
200 
201  return true;
202 }
203 
204 
205 // ************************************************************************* //
scalar delta
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
virtual bool write()
Calculate the histogram and write.
Definition: histogram.C:118
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
static bool master(const label communicator=0)
Am I the master process.
Definition: UPstream.H:411
word format(conversionProperties.lookup("format"))
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:715
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:68
virtual bool execute()
Execute, currently does nothing.
Definition: histogram.C:112
Macros for easy insertion into run-time selection tables.
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
static List< word > fieldNames
Definition: globalFoam.H:46
bool read(const char *, int32_t &)
Definition: int32IO.C:85
virtual bool read(const dictionary &)
Read the histogram data.
Definition: histogram.C:98
dynamicFvMesh & mesh
Holds list of sampling positions.
Definition: coordSet.H:49
A class for handling words, derived from string.
Definition: word.H:59
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
static autoPtr< writer > New(const word &writeFormat)
Return a reference to the selected writer.
Definition: writer.C:35
static const char nl
Definition: Ostream.H:262
bool mkDir(const fileName &, mode_t=0777)
Make a directory and return an error if it could not be created.
Definition: POSIX.C:295
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
List< word > wordList
A List of words.
Definition: fileName.H:54
vector point
Point is a vector.
Definition: point.H:41
fileName::Type type(const fileName &)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:461
addToRunTimeSelectionTable(functionObject, blendingFactor, dictionary)
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
const DimensionedField< scalar, volMesh > & V() const
Return cell volumes.
#define Log
Report write to Foam::Info if the local log switch is true.
dimensioned< scalar > mag(const dimensioned< Type > &)
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:53
defineTypeNameAndDebug(fvMeshFunctionObject, 0)
functionObject base class for writing single files
Definition: writeFile.H:56
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
const word & name() const
Return name.
Definition: IOobject.H:260
Namespace for OpenFOAM.
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:243
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:451