Go to the documentation of this file.
35 namespace distributionModels
49 void Foam::distributionModels::binned::initialise()
51 const label nSample(xy_.size());
54 for (label bini = 1; bini < nSample; ++bini)
56 xy_[bini][1] += xy_[bini - 1][1];
60 scalar sumProb = xy_.last()[1];
65 <<
type() <<
"distribution: "
66 <<
"The sum of elements in the second column cannot be zero." <<
nl
67 <<
"sum = " << sumProb
73 xy_[bini][1] /= sumProb;
87 meanValue_ = xy_[bini][1];
100 xy_(distributionModelDict_.lookup(
"distribution")),
103 minValue_ = xy_[0][0];
104 maxValue_ = xy_[xy_.size()-1][0];
115 const scalar binWidth,
127 minValue_ =
min(minValue_, sampleData[i]);
128 maxValue_ =
max(maxValue_, sampleData[i]);
131 const label bin0 = floor(minValue_/binWidth);
132 const label bin1 = ceil(maxValue_/binWidth);
133 const label nBin = bin1 - bin0;
138 <<
"Data cannot be binned - zero bins generated" <<
nl
139 <<
" Bin width : " << binWidth <<
nl
140 <<
" Sample data : " << sampleData
150 xy_[bini][0] = (bin0 + bini)*binWidth;
158 label bini = floor(sampleData[i]/binWidth) - bin0;
159 label binii =
min(bini + 1, nBin - 1);
161 scalar d1 =
mag(sampleData[i] - xy_[bini][0]);
162 scalar d2 =
mag(xy_[binii][0] - sampleData[i]);
182 meanValue_(
p.meanValue_)
190 const scalar u = rndGen_.sample01<scalar>();
192 for (label i = 0; i < xy_.size() - 1; ++i)
231 dict.add(
"distribution", xy_);
240 dict.readEntry(
"distribution", xy_);
virtual void readDict(const dictionary &dict)
Read data from dictionary.
virtual void writeData(Ostream &os) const
Write data to stream.
A class for handling words, derived from Foam::string.
const word dictName("faMeshDefinition")
virtual void readData(Istream &os)
Read data from stream.
addToRunTimeSelectionTable(distributionModel, binned, dictionary)
Particle-size distribution model wherein random samples are drawn from a given discrete set of (bin,...
Istream & operator>>(Istream &, directionInfo &)
Ostream & endl(Ostream &os)
Add newline and flush stream.
A library of runtime-selectable doubly-truncated probability distribution models. Returns random samp...
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.
static const dictionary null
An empty dictionary, which is also the parent for all dictionaries.
virtual scalar meanValue() const
Return the arithmetic mean of the distribution data.
virtual scalar sample() const
Sample the distribution.
static void check(const int retVal, const char *what)
const dimensionedScalar b
Wien displacement law constant: default SI units: [m.K].
An Istream is an abstract base class for all input systems (streams, files, token lists etc)....
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
virtual bool check(const char *operation) const
Check IOstream status for given operation.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
OBJstream os(runTime.globalPath()/outputName)
Macros for easy insertion into run-time selection tables.
static const char * header
errorManipArg< error, int > exit(error &err, const int errNo=1)
binned(const dictionary &dict, Random &rndGen)
Construct from dictionary.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
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.
defineTypeNameAndDebug(binned, 0)
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
virtual dictionary writeDict(const word &dictName) const
Write data in dictionary format.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
#define WarningInFunction
Report a warning using Foam::Warning.