Go to the documentation of this file.
36 namespace distributionModels
46 void Foam::distributionModels::general::initialise()
48 static scalar eps = ROOTVSMALL;
57 for (label i = 1; i < nEntries_; ++i)
59 scalar
k = (xy_[i][1] - xy_[i-1][1])/(xy_[i][0] - xy_[i-1][0] + eps);
60 scalar d = xy_[i-1][1] -
k*xy_[i-1][0];
61 scalar
y1 = xy_[i][0]*(0.5*
k*xy_[i][0] + d);
62 scalar
y0 = xy_[i-1][0]*(0.5*
k*xy_[i-1][0] + d);
67 integral_[i] = xy_[i][1];
72 integral_[i] =
area + integral_[i-1];
74 y1 =
sqr(xy_[i][0])*(1.0/3.0*
k*xy_[i][0] + 0.5*d);
75 y0 =
sqr(xy_[i-1][0])*(1.0/3.0*
k*xy_[i-1][0] + 0.5*d);
76 meanValue_ +=
y1 -
y0;
81 const scalar sumArea = integral_.last();
83 for (label i = 0; i < nEntries_; ++i)
85 xy_[i][1] /= sumArea + eps;
86 integral_[i] /= sumArea + eps;
89 meanValue_ /= sumArea + eps;
90 meanValue_ = cumulative_ ? (
maxValue_ - meanValue_) : meanValue_;
103 xy_(distributionModelDict_.lookup(
"distribution")),
104 nEntries_(xy_.size()),
106 integral_(nEntries_),
107 cumulative_(distributionModelDict_.getOrDefault(
"cumulative",
false))
109 minValue_ = xy_[0][0];
110 maxValue_ = xy_[nEntries_-1][0];
115 if (cumulative_ && xy_[0][1] != 0)
118 <<
type() <<
"distribution: "
119 <<
"Elements in the second column for cumulative "
120 <<
"distribution functions must start from zero." <<
nl
121 <<
"First element = " << xy_[0][1]
125 for (label i = 0; i < nEntries_; ++i)
127 if (i > 0 && xy_[i][0] <= xy_[i-1][0])
130 <<
type() <<
"distribution: "
131 <<
"Elements in the first column must "
132 <<
"be specified in an ascending order." <<
nl
133 <<
"Please see the row i = " << i <<
nl
134 <<
"xy[i-1] = " << xy_[i-1] <<
nl
135 <<
"xy[i] = " << xy_[i]
139 if (xy_[i][0] < 0 || xy_[i][1] < 0)
142 <<
type() <<
"distribution: "
143 <<
"Input pairs cannot contain any negative element." <<
nl
144 <<
"Please see the row i = " << i <<
nl
145 <<
"xy[i] = " << xy_[i]
149 if (cumulative_ && i > 0 && xy_[i][1] < xy_[i-1][1])
152 <<
type() <<
"distribution: "
153 <<
"Elements in the second column for cumulative "
154 <<
"distribution functions must be non-decreasing." <<
nl
155 <<
"Please see the row i = " << i <<
nl
156 <<
"xy[i-1] = " << xy_[i-1] <<
nl
157 <<
"xy[i] = " << xy_[i]
169 const scalar binWidth,
184 minValue_ =
min(minValue_, sampleData[i]);
185 maxValue_ =
max(maxValue_, sampleData[i]);
188 label bin0 = floor(minValue_/binWidth);
189 label bin1 = ceil(maxValue_/binWidth);
190 nEntries_ = bin1 - bin0;
195 <<
"Data cannot be binned - zero bins generated" <<
nl
196 <<
" Bin width : " << binWidth <<
nl
197 <<
" Sample data : " << sampleData
203 xy_.setSize(nEntries_);
206 for (label bini = 0; bini < nEntries_; ++bini)
208 xy_[bini][0] = (bin0 + bini)*binWidth;
215 label bini = floor(sampleData[i]/binWidth) - bin0;
227 nEntries_(
p.nEntries_),
228 meanValue_(
p.meanValue_),
229 integral_(
p.integral_),
230 cumulative_(
p.cumulative_)
238 const scalar u = rndGen_.sample01<scalar>();
242 while (integral_[
n] <= u)
247 const scalar
k = (xy_[
n][1] - xy_[
n-1][1])/(xy_[
n][0] - xy_[
n-1][0]);
248 const scalar d = xy_[
n-1][1] -
k*xy_[
n-1][0];
256 u + xy_[
n-1][0]*(0.5*
k*xy_[
n-1][0] + d) - integral_[
n-1];
261 const scalar
p = 2.0*d/
k;
262 const scalar q = -2.0*
alpha/
k;
263 const scalar sqrtEr =
sqrt(0.25*
p*
p - q);
265 const scalar x1 = -0.5*
p + sqrtEr;
266 const scalar x2 = -0.5*
p - sqrtEr;
267 if ((x1 >= xy_[
n-1][0]) && (x1 <= xy_[
n][0]))
322 xy_.setSize(
x.size());
A class for handling words, derived from Foam::string.
virtual tmp< scalarField > x() const
Bin boundaries.
dimensionedScalar y1(const dimensionedScalar &ds)
virtual scalar meanValue() const
Return the arithmetic mean of the distribution data.
A class for managing temporary objects.
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
const word dictName("faMeshDefinition")
addToRunTimeSelectionTable(distributionModel, binned, dictionary)
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.
scalar maxValue_
Maximum of the distribution.
Particle-size distribution model wherein random samples are drawn from a given arbitrary probability ...
#define forAll(list, i)
Loop across all elements in list.
virtual tmp< scalarField > y() const
Probabilities.
static const dictionary null
An empty dictionary, which is also the parent for all dictionaries.
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].
General relative velocity model.
An Istream is an abstract base class for all input systems (streams, files, token lists etc)....
void setSize(const label n)
Alias for resize()
virtual void readDict(const dictionary &dict)
Read data from dictionary.
dimensionedScalar y0(const dimensionedScalar &ds)
virtual dictionary writeDict(const word &dictName) const
Write data in dictionary format.
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.
general(const dictionary &dict, Random &rndGen)
Construct from components.
errorManipArg< error, int > exit(error &err, const int errNo=1)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensionedScalar sqrt(const dimensionedScalar &ds)
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...
label k
Boltzmann constant.
virtual void readData(Istream &os)
Read data from stream.
static tmp< T > New(Args &&... args)
Construct tmp of T with forwarding arguments.
const wordList area
Standard area field types (scalar, vector, tensor, etc)
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
virtual void writeData(Ostream &os) const
Write data to stream.
#define WarningInFunction
Report a warning using Foam::Warning.