Go to the documentation of this file.
37 namespace functionObjects
47 void Foam::functionObjects::histogram::writeGraph
49 const coordSet& coords,
50 const word& fieldName,
60 /formatterPtr_().getFileName
67 Log <<
" Writing histogram of " << fieldName
68 <<
" to " << graphFile.name() <<
endl;
71 fieldNames[0] = fieldName;
72 fieldNames[1] = fieldName +
"Count";
73 List<const scalarField*> yPtrs(2);
74 yPtrs[0] = &normalizedValues;
75 yPtrs[1] = &absoluteValues;
76 formatterPtr_().write(coords, fieldNames, yPtrs, graphFile);
105 dict.readEntry(
"field", fieldName_);
107 max_ =
dict.getOrDefault<scalar>(
"max", -GREAT);
108 min_ =
dict.getOrDefault<scalar>(
"min", GREAT);
109 dict.readEntry(
"nBins", nBins_);
114 <<
"Number of histogram bins = " << nBins_
115 <<
" cannot be negative or zero."
139 Log <<
" Looking up field " << fieldName_ <<
endl;
143 Log <<
" Reading field " << fieldName_ <<
endl;
151 mesh_.time().timeName(),
169 scalar histMax = max_;
170 scalar histMin = min_;
181 Log <<
" Determined histogram bounds from field"
182 <<
" min/max(" << fieldName_ <<
") = "
183 << histMin <<
' ' << histMax <<
endl;
185 else if (min_ == GREAT)
192 const scalar
delta = (histMax- histMin)/nBins_;
194 scalar
x = histMin + 0.5*
delta;
207 const label bini = (
field[celli] - histMin)/
delta;
208 if (bini >= 0 && bini < nBins_)
210 dataNormalized[bini] += V[celli];
220 const scalar sumData =
sum(dataNormalized);
224 dataNormalized /= sumData;
239 count[i] = 1.0*dataCount[i];
242 writeGraph(coords, fieldName_, dataNormalized,
count);
void reset(T *p=nullptr) noexcept
Delete managed object and set to new given pointer.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
A class for handling words, derived from Foam::string.
static constexpr const zero Zero
Global zero (0)
virtual bool write()
Calculate the histogram and write.
bool read(const char *buf, int32_t &val)
Same as readInt32.
bool valid() const noexcept
True if the managed pointer is non-null.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Specialization of Foam::functionObject for an Foam::fvMesh, providing a reference to the Foam::fvMesh...
virtual bool execute()
Execute (effectively no-op)
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
#define forAll(list, i)
Loop across all elements in list.
List< word > wordList
A List of words.
word format(conversionProperties.get< word >("format"))
virtual bool read(const dictionary &dict)
Read.
word name(const complex &c)
Return string representation of complex.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
virtual bool read(const dictionary &dict)
Read optional controls.
Macros for easy insertion into run-time selection tables.
Holds list of sampling positions.
errorManip< error > abort(error &err)
virtual bool read(const dictionary &)
Read the histogram data.
static bool master(const label communicator=0)
Am I the master process.
static void listCombineGather(const List< commsStruct > &comms, List< T > &Value, const CombineOp &cop, const int tag, const label comm)
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
addToRunTimeSelectionTable(functionObject, ObukhovLength, dictionary)
unsigned int count(const UList< bool > &bools, const bool val=true)
Count number of 'true' entries.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
defineTypeNameAndDebug(ObukhovLength, 0)
histogram(const word &name, const Time &runTime, const dictionary &dict)
Construct from Time and dictionary.
Base class for writing single files from the function objects.
vector point
Point is a vector.
fileName baseTimeDir() const
Return the base directory for the current time value.
bool mkDir(const fileName &pathName, mode_t mode=0777)
Make a directory and return an error if it could not be created.
static autoPtr< writer > New(const word &writeFormat)
Return a reference to the selected writer.