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 | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2016 OpenFOAM Foundation
9 Copyright (C) 2016-2022 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27\*---------------------------------------------------------------------------*/
28
29#include "histogram.H"
30#include "volFields.H"
31#include "ListOps.H"
33
34// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
35
36namespace Foam
37{
38namespace functionObjects
39{
42}
43}
44
45
46// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
47
49(
50 const word& name,
51 const Time& runTime,
52 const dictionary& dict
53)
54:
55 functionObjects::fvMeshFunctionObject(name, runTime, dict),
56 functionObjects::writeFile(obr_, name),
57 max_(-GREAT),
58 min_(GREAT),
59 setWriterPtr_(nullptr)
60{
61 read(dict);
62}
63
64
65// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
66
68{
71
72 dict.readEntry("field", fieldName_);
73
74 max_ = dict.getOrDefault<scalar>("max", -GREAT);
75 min_ = dict.getOrDefault<scalar>("min", GREAT);
76 dict.readEntry("nBins", nBins_);
77
78 if (nBins_ < 1)
79 {
81 << "Number of histogram bins = " << nBins_
82 << " cannot be negative or zero."
83 << abort(FatalError);
84 }
85
86 const word writeType(dict.get<word>("setFormat"));
87
88 setWriterPtr_ = coordSetWriter::New
89 (
90 writeType,
91 dict.subOrEmptyDict("formatOptions").optionalSubDict(writeType)
92 );
93
94 return true;
95}
96
97
99{
100 return true;
101}
102
103
105{
106 Log << type() << " " << name() << " write:" << nl;
107
108 tmp<volScalarField> tfield;
109 tfield.cref(obr_.cfindObject<volScalarField>(fieldName_));
110
111 if (tfield)
112 {
113 Log << " Looking up field " << fieldName_ << endl;
114 }
115 else
116 {
117 Log << " Reading field " << fieldName_ << endl;
119 (
121 (
122 fieldName_,
123 mesh_.time().timeName(),
124 mesh_,
127 ),
128 mesh_
129 );
130 }
131 const auto& field = tfield();
132
133 scalar histMax = max_;
134 scalar histMin = min_;
135
136 if (max_ == -GREAT)
137 {
138 // Determine current min and max
139 histMax = max(field).value();
140
141 if (min_ == GREAT)
142 {
143 histMin = min(field).value();
144 }
145 Log << " Determined histogram bounds from field"
146 << " min/max(" << fieldName_ << ") = "
147 << histMin << ' ' << histMax << endl;
148 }
149 else if (min_ == GREAT)
150 {
151 histMin = 0;
152 }
153
154 // Calculate the mid-points of bins for the graph axis
155 pointField xBin(nBins_, Zero);
156 const scalar delta = (histMax - histMin)/nBins_;
157
158 {
159 scalar x = histMin + 0.5*delta;
160 for (point& p : xBin)
161 {
162 p.x() = x;
163 x += delta;
164 }
165 }
166
167 scalarField dataNormalized(nBins_, Zero);
168 labelField dataCount(nBins_, Zero);
169 const scalarField& V = mesh_.V();
170
171 forAll(field, celli)
172 {
173 const label bini = (field[celli] - histMin)/delta;
174 if (bini >= 0 && bini < nBins_)
175 {
176 dataNormalized[bini] += V[celli];
177 dataCount[bini]++;
178 }
179 }
180
183
184 if (Pstream::master())
185 {
186 const scalar sumData = sum(dataNormalized);
187
188 if (sumData > SMALL)
189 {
190 dataNormalized /= sumData;
191
192 const coordSet coords(fieldName_, "x", xBin, mag(xBin));
193
194 auto& writer = *setWriterPtr_;
195
196 writer.open
197 (
198 coords,
199 (
201 / (coords.name() + coordSetWriter::suffix(fieldName_))
202 )
203 );
204
205 Log << " Writing histogram of " << fieldName_
206 << " to " << writer.path() << endl;
207
208 writer.nFields(2);
209 writer.write(fieldName_, dataNormalized);
210 writer.write(fieldName_ + "Count", dataCount);
211
212 writer.close(true);
213 }
214 }
215
216 return true;
217}
218
219
220// ************************************************************************* //
scalar delta
Various functions to operate on Lists.
#define Log
Definition: PDRblock.C:35
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
vtk::internalMeshWriter writer(topoMesh, topoCells, vtk::formatType::INLINE_ASCII, runTime.path()/"blockTopology")
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:170
static void listCombineGather(const List< commsStruct > &comms, List< T > &values, const CombineOp &cop, const int tag, const label comm)
virtual bool read()
Re-read model coefficients if they have changed.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:80
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
Definition: Time.C:717
static word suffix(const word &fldName, const word &fileExt=word::null)
Name suffix based on fieldName (underscore separator)
Holds list of sampling positions.
Definition: coordSet.H:56
const word & name() const noexcept
The coord-set name.
Definition: coordSet.H:134
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
Abstract base-class for Time/database function objects.
Watches for presence of the named trigger file in the case directory and signals a simulation stop (o...
Definition: abort.H:128
Specialization of Foam::functionObject for an Foam::fvMesh, providing a reference to the Foam::fvMesh...
Computes the volume-weighted histogram of an input volScalarField.
Definition: histogram.H:188
virtual bool execute()
Execute (effectively no-op)
Definition: histogram.C:98
virtual bool write()
Calculate the histogram and write.
Definition: histogram.C:104
virtual bool read(const dictionary &)
Read the histogram data.
Definition: histogram.C:67
Computes the magnitude of an input field.
Definition: mag.H:153
Base class for writing single files from the function objects.
Definition: writeFile.H:120
fileName baseTimeDir() const
Return the base directory for the current time value.
Definition: writeFile.C:73
splitCell * master() const
Definition: splitCell.H:113
A class for managing temporary objects.
Definition: tmp.H:65
const T & cref() const
Definition: tmpI.H:213
A class for handling words, derived from Foam::string.
Definition: word.H:68
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
Definition: className.H:121
volScalarField & p
rDeltaTY field()
engineTime & runTime
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
Namespace for OpenFOAM.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:47
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
Definition: MSwindows.C:598
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:33
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
error FatalError
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:59
constexpr char nl
The newline '\n' character (0x0a)
Definition: Ostream.H:53
dictionary dict
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333