 |
The SKIRT project
advanced radiative transfer for astrophysics
|
|
The SKIRT project
advanced radiative transfer for astrophysics
|
#include <Random.hpp>
Public Member Functions | |
| double | cdfLinLin (const Array &xv, const Array &Pv) |
| double | cdfLogLog (const Array &xv, const Array &pv, const Array &Pv) |
| Direction | direction () |
| Direction | direction (Direction bfk, double costheta) |
| double | expon () |
| double | exponCutoff (double xmax) |
| double | gauss () |
| double | lorentz () |
| Vec | maxwell () |
| void | pop () |
| Position | position (const Box &box) |
| void | push (int seed) |
| int | seed () const |
| double | uniform () |
 Public Member Functions inherited from SimulationItem | |
| template<class T > | |
| T * | find (bool setup=true) const |
| template<class T > | |
| T * | interface (int levels=-999999, bool setup=true) const |
| virtual string | itemName () const |
| void | setup () |
| string | typeAndName () const |
| Public Member Functions inherited from Item | |
| Item (const Item &)=delete | |
| virtual | ~Item () |
| void | addChild (Item *child) |
| const vector< Item * > & | children () const |
| virtual void | clearItemListProperty (const PropertyDef *property) |
| void | destroyChild (Item *child) |
| virtual bool | getBoolProperty (const PropertyDef *property) const |
| virtual vector< double > | getDoubleListProperty (const PropertyDef *property) const |
| virtual double | getDoubleProperty (const PropertyDef *property) const |
| virtual string | getEnumProperty (const PropertyDef *property) const |
| virtual int | getIntProperty (const PropertyDef *property) const |
| virtual vector< Item * > | getItemListProperty (const PropertyDef *property) const |
| virtual Item * | getItemProperty (const PropertyDef *property) const |
| virtual string | getStringProperty (const PropertyDef *property) const |
| int | getUtilityProperty (string name) const |
| virtual void | insertIntoItemListProperty (const PropertyDef *property, int index, Item *item) |
| Item & | operator= (const Item &)=delete |
| Item * | parent () const |
| virtual void | removeFromItemListProperty (const PropertyDef *property, int index) |
| virtual void | setBoolProperty (const PropertyDef *property, bool value) |
| virtual void | setDoubleListProperty (const PropertyDef *property, vector< double > value) |
| virtual void | setDoubleProperty (const PropertyDef *property, double value) |
| virtual void | setEnumProperty (const PropertyDef *property, string value) |
| virtual void | setIntProperty (const PropertyDef *property, int value) |
| virtual void | setItemProperty (const PropertyDef *property, Item *item) |
| virtual void | setStringProperty (const PropertyDef *property, string value) |
| void | setUtilityProperty (string name, int value) |
| virtual string | type () const |
Protected Member Functions | |
| Random () | |
| void | setupSelfBefore () override |
| Protected Member Functions inherited from SimulationItem | |
| SimulationItem () | |
| virtual bool | offersInterface (const std::type_info &interfaceTypeInfo) const |
| virtual void | setupSelfAfter () |
| virtual void | setupSelfBefore () |
| Protected Member Functions inherited from Item | |
| Item () | |
Private Types | |
| using | BaseType = SimulationItem |
| using | ItemType = Random |
Private Attributes | |
| int | _seed |
Friends | |
| class | ItemRegistry |
The Random class offers pseudo-random number generation facilities in a multi-threaded environment. At the core of its operation is the function to generate a uniform deviate, i.e. a pseudo-random number drawn from the uniform probability distribution over the unit interval. Additional functions allow drawing pseudo-random numbers from other frequently-used probability distributions.
To avoid the need for synchronization between multiple execution threads, each thread receives its own thread-local pseudo-random number generator instance. The parent thread in each process (i.e. the thread that calls the setup() function) receives a predictable generator. This generator is initialized with a fixed state depending only on the value of the user-configurable seed property, so that it delivers exactly the same pseudo-random sequence for every parent thread in every process and for every execution of the program. All other (child) threads receive an arbitrary generator. These generators are seeded with a truly random state obtained from the operating system, so that they deliver a unique and unpredictable pseudo-random sequence for each thread, and for each execution of the program.
This mechanism is intended to support the following use cases. In single thread/single process mode, repeated execution (with the same seed value) must produce identical results. This is important, for example, for automated testing. The user can force a different randomization by specifying another seed value. In multi-threading and/or multi-processing mode, for tasks employing parallelization, each of the threads in each of the processes requires a different random sequence. On the other hand, in this mode, for serialized tasks, the (single) thread employed in each of the processes must receive the same random sequence. This is important for functions that rely on randomness but need to produce the same result in every parallel process.
The recommended use of the Random class is to include a single instance in each simulation run-time hierarchy. The use cases discussed above can be achieved as follows. In single thread/single process mode, perform all tasks in the parent thread. In multi-threading and/or multi-processing mode, perform serial tasks in the parent thread of each process, and perform all parallized tasks in a child thread.
As an additional service, this class allows temporarily installing a predictable random number generator with a given seed in the current execution thread through the push() and pop() functions. This supports the use case where, usually during setup, the same pseudo-random sequence is required in multiple places.
All random number generators used in this class are based on the 64-bit Mersenne twister, which offers a sufficiently long period and acceptable spectral properties for most purposes.
|
inlineprotected |
This function generates a random number drawn from an arbitrary probability distribution \(p(x)\,{\text{d}}x\) with corresponding cumulative distribution function \(P(x)\). The function accepts a discretized version \(P_i\) of the cdf sampled at a set of \(N\) points \(x_i\). A uniform deviate \(\cal{X}\) is generated, and the equation \({\cal{X}}=P(x)\) is solved using linear interpolation (i.e. assuming piece-wise linear behavior of the cdf (and equivalently, of the underlying pdf).
This function generates a random number drawn from an arbitrary probability distribution \(p(x)\,{\text{d}}x\) with corresponding cumulative distribution function \(P(x)\). The function accepts discretized versions \(p_i\) and \(P_i\) of the pdf and cdf sampled at a set of \(N\) points \(x_i\). A uniform deviate \(\cal{X}\) is generated, and the equation \({\cal{X}}=P(x)\) is solved using a interpolation mechanism that assumes that the pdf is linear in log-log space between any two grid points (equivalent to power-law behavior), as described below.
Consider the pdf values \(p_i\) and \(p_{i+1}\) at two consecutive grid points \(x_i\) and \(x_{i+1}\). Assuming power-law behavior, the pdf between these two grid points can be written as
\[ p(x) = p_i \left(\frac{x}{x_i}\right)^{\alpha_i}, \quad\mathrm{with}\; \alpha_i = \frac{\ln p_{i+1}/\ln p_i}{\ln x_{i+1}/\ln x_i} \]
With \(\mathcal{X}\) a random deviate for which \(x_i\leq\mathcal{X}\leq x_{i+1}\) happens to be true, we thus need to invert the relation
\[ \mathcal{X}-x_i = \int_{x_i}^x p(x')\,\mathrm{d}x' = \int_{x_i}^x p_i \left(\frac{x'}{x_i}\right)^{\alpha_i}\mathrm{d}x' = p_i x_i \;\mathrm{gln}\left(-\alpha_i, \frac{x}{x_i}\right) \]
which leads to
\[x = x_i \;\mathrm{gexp}\left(-\alpha_i, \frac{\mathcal{X}-x_i}{p_i x_i}\right)\]
where \(\mathrm{gln}(a,x)\) is the generalized logarithm and \(\mathrm{gexp}(a,x)\) the generalized exponential, defined in the description of respectively the SpecialFunctions::gln() and SpecialFunctions::gexp() functions.
| Direction Random::direction | ( | ) |
This function generates a random direction on the unit sphere, i.e. a couple \((\theta,\phi)\) from the two-dimensional probability density
\[ p(\theta,\phi)\,d\theta\,d\phi = \left(\frac{\sin\theta}{2}\,d\theta\right) \left(\frac{1}{2\pi}\,d\varphi\right).\]
A random direction on the unit sphere can thus be constructed by taking two uniform deviates \({\cal{X}}_1\) and \({\cal{X}}_2\), and solving the two equations
\[ \begin{split} {\cal{X}}_1 &= \int_0^\theta \frac{\sin\theta'\,d\theta'}{2} \\ {\cal{X}}_2 &= \int_0^\varphi \frac{d\varphi'}{2\pi} \end{split} \]
for \(\theta\) and \(\varphi\). The solution is readily found,
\[ \begin{split} \theta &= \arccos\left(2{\cal{X}}_1-1\right) \\ \varphi &= 2\pi\,{\cal {X}}_2. \end{split} \]
Once these spherical coordinates are calculated, a Direction object can be constructed by calling the constructor Direction::Direction(double theta, double phi).
This function generates a new direction on the unit sphere deviating from a given original direction \(\bf{k}\) by a given polar angle \(\theta\) (specified through its cosine) and a uniformly random azimuth angle \(\phi\). The function can use an arbitrary reference point for the azimuth angle since it is distributed uniformly. We use equation (A37) of Bianchi et al. 1996 (ApJ 465,127), which uses the projection of the z-axis as the reference point.
| double Random::expon | ( | ) |
This function generates a random number from an exponential distribution function, defined by the probability distribution
\[ p(x)\,{\rm d}x = {\rm e}^{-x}\,{\rm d}x.\]
A simple inversion technique is used.
| double Random::exponCutoff | ( | double | xmax | ) |
This function generates a random number from an exponential distribution function with a cut-off, defined by the probability distribution
\[ p(x)\,{\rm d}x \propto \begin{cases} \;{\rm e}^{-x}\,{\rm d}x &\qquad \text{for $0<x<x_{\text{max}}$,} \\ \;0 &\qquad \text{for $x>x_{\text{max}}$.} \end{cases} \]
A simple inversion technique is used.
| double Random::gauss | ( | ) |
This function generates a random number from a Gaussian distribution function with mean 0 and standard deviation 1, i.e. defined by the probability distribution
\[ p(x)\,{\rm d}x = \frac{1}{\sqrt{2\pi}}\, {\rm e}^{-\frac12\,x^2}\,{\rm d}x.\]
The algorithm used and the implementation are taken from Press et al. (2002).
| double Random::lorentz | ( | ) |
This function generates a random number from a Lorentzian (or Cauchy) distribution function with its peak at 0 and a half-width at half-maximum (HWHM) of 1, i.e. defined by the probability distribution
\[ p(x)\,{\rm d}x = \frac{1}{\pi}\, \frac{1}{x^2+1}\,{\rm d}x.\]
The cumulative distribution function is easily calculated and inverted analytically, so that a simple inversion techique can be used.
| Vec Random::maxwell | ( | ) |
This function generates a random velocity from a three-dimensional Maxwell-Boltzmann distribution with velocity dispersion 1, which is equivalent to a Gaussian distribution with mean 0 and standard deviation 1 for each of the Cartesian velocity components. The implementation simply calls the function gauss() to obtain each of these three velocity components.
| void Random::pop | ( | ) |
This function pops the most recently pushed random number generator from the stack for the current thread and establishes it as the active random number generator for the current thread. If the stack does not contain a random number generator, the behavior of this function is undefined.
This function generates a uniformly distributed random position in a given box (i.e. a cuboid lined up with the coordinate axes).
| void Random::push | ( | int | seed | ) |
This function pushes the active random number generator for the current thread onto the stack for the current thread and establishes a new predictable random number generator, initialized with the specified seed, as the active random number generator for the current thread.
|
inline |
This function returns the value of the discoverable integer property seed : "the seed for the random generator" .
The minimum value for this property is "0" .
The maximum value for this property is "1000000" .
The default value for this property is given by the conditional value expression "0" .
This property is displayed only if the Boolean expression "Level3" evaluates to true after replacing the names by true or false depending on their presence.
|
overrideprotectedvirtual |
This function initializes the random generator for the calling thread (deemed the parent thread) to its fixed initial state depending on the value of the user-configurable seed property.
Reimplemented from SimulationItem.
| double Random::uniform | ( | ) |
This function generates a uniform deviate, i.e. a random double precision number in the open interval (0,1). The interval borders zero and one are never returned.