diff --git a/Repository/RandomGenerator.cc b/Repository/RandomGenerator.cc --- a/Repository/RandomGenerator.cc +++ b/Repository/RandomGenerator.cc @@ -1,174 +1,184 @@ // -*- C++ -*- // // RandomGenerator.cc is a part of ThePEG - Toolkit for HEP Event Generation // Copyright (C) 1999-2017 Leif Lonnblad // // ThePEG is licenced under version 3 of the GPL, see COPYING for details. // Please respect the MCnet academic guidelines, see GUIDELINES for details. // // // This is the implementation of the non-inlined, non-templated member // functions of the RandomGenerator class. // #include "RandomGenerator.h" #include "ThePEG/Persistency/PersistentOStream.h" #include "ThePEG/Persistency/PersistentIStream.h" #include "ThePEG/Interface/Parameter.h" #include "ThePEG/Interface/ClassDocumentation.h" #include "gsl/gsl_randist.h" extern "C" { typedef struct { ThePEG::RandomGenerator * r; } thepeg_random_state_t; void thepeg_random_set(void *, unsigned long int) {} double thepeg_random_get_double(void * s) { return static_cast(s)->r->rnd(); } unsigned long int thepeg_random_get(void * s) { return static_cast(std::numeric_limits::max()*thepeg_random_get_double(s)); } static const gsl_rng_type thepeg_random_type = {"thepeg_random", (unsigned long int)std::numeric_limits::max(), 0, sizeof(thepeg_random_state_t), &thepeg_random_set, &thepeg_random_get, &thepeg_random_get_double}; const gsl_rng_type *gsl_rng_thepeg_random = &thepeg_random_type; } using namespace ThePEG; RandomGenerator::RandomGenerator() : theNumbers(1000), theSize(1000), theSeed(-1), savedGauss(0.0), gaussSaved(false) { nextNumber = theNumbers.end(); gsl = gsl_rng_alloc(gsl_rng_thepeg_random); static_cast(gsl->state)->r = this; } RandomGenerator::RandomGenerator(const RandomGenerator & rg) : Interfaced(rg), theNumbers(rg.theNumbers), theSize(rg.theSize), theSeed(rg.theSeed), savedGauss(rg.savedGauss), gaussSaved(rg.gaussSaved) { nextNumber = theNumbers.begin() + ( RndVector::const_iterator(rg.nextNumber) - rg.theNumbers.begin() ); gsl = gsl_rng_alloc(gsl_rng_thepeg_random); static_cast(gsl->state)->r = this; } RandomGenerator::~RandomGenerator() { gsl_rng_free(gsl); } void RandomGenerator::doinit() { if ( theSeed != 0 ) setSeed(theSeed); flush(); } void RandomGenerator::setSize(size_type newSize) { RndVector newNumbers(newSize); RndVector::iterator nextNew = newNumbers.end() - min( int(theNumbers.end() - nextNumber), int(newSize) ); for ( RndVector::iterator i = nextNew; i != newNumbers.end(); ++i ) *i = *nextNumber++; RndVector::difference_type pos = nextNew - newNumbers.begin(); theNumbers.swap(newNumbers); nextNumber = theNumbers.begin() + pos; } bool RandomGenerator::prndbool(double p) { if ( p >= 1.0 ) return true; if ( p <= 0.0 ) return false; double r = rnd(); if ( r < p ) { push_back(r/p); return true; } else { push_back((r - p)/(1.0 - p)); return false; } } int RandomGenerator::rndsign(double p1, double p2, double p3) { double sum = p1 + p2 + p3; double r = rnd()*sum; if ( r < p1 ) return -1; else if ( r < p1 + p2 ) return 0; else return 1; } int RandomGenerator::prndsign(double p1, double p2, double p3) { double sum = p1 + p2 + p3; double r = rnd()*sum; if ( r < p1 ) { push_back(r/p1); return -1; } else if ( r < p1 + p2 ) { push_back((r - p1)/p2); return 0; } else { push_back((r - p1 - p2)/p3); return 1; } } int RandomGenerator::rnd4(double p0, double p1, double p2, double p3) { double sum = p0 + p1 + p2 + p3; double r = rnd()*sum; if ( r < p0 ) return 0; else if ( r < p0 + p1 ) return 1; else if ( r < p0 + p1 + p2 ) return 2; else return 3; } +int RandomGenerator::rnd5(double p0, double p1, double p2, double p3, double p4) { + double sum = p0 + p1 + p2 + p3 + p4; + double r = rnd()*sum; + if ( r < p0 ) return 0; + else if ( r < p0 + p1 ) return 1; + else if ( r < p0 + p1 + p2 ) return 2; + else if ( r < p0 + p1 + p2 + p3 ) return 3; + else return 4; +} + long RandomGenerator::rndPoisson(double mean) { return gsl_ran_poisson(gsl, mean); } void RandomGenerator::persistentOutput(PersistentOStream & os) const { os << theNumbers << RndVector::const_iterator(nextNumber) - theNumbers.begin() << theSize << theSeed << savedGauss << gaussSaved; } void RandomGenerator::persistentInput(PersistentIStream & is, int) { RndVector::difference_type pos; is >> theNumbers >> pos >> theSize >> theSeed >> savedGauss >> gaussSaved; nextNumber = theNumbers.begin() + pos; } ClassDescription RandomGenerator::initRandomGenerator; void RandomGenerator::Init() { static ClassDocumentation documentation ("There is no documentation for the ThePEG::RandomGenerator class"); static Parameter interfaceSize ("CacheSize", "The Random numbers are generated in chunks of this size.", &RandomGenerator::theSize, 1000, 10, 100000, true, false, true, &RandomGenerator::setSize); static Parameter interfaceSeed ("Seed", "The seed with which this random generator is initialized. " "If set to -1, the default build-in seed will be used. If set to zero, no seed will " "be set.", &RandomGenerator::theSeed, -1, -1, 100000000, true, false, false); interfaceSeed.setHasDefault(false); interfaceSize.rank(10); interfaceSeed.rank(9); } diff --git a/Repository/RandomGenerator.h b/Repository/RandomGenerator.h --- a/Repository/RandomGenerator.h +++ b/Repository/RandomGenerator.h @@ -1,491 +1,497 @@ // -*- C++ -*- // // RandomGenerator.h is a part of ThePEG - Toolkit for HEP Event Generation // Copyright (C) 1999-2017 Leif Lonnblad // // ThePEG is licenced under version 3 of the GPL, see COPYING for details. // Please respect the MCnet academic guidelines, see GUIDELINES for details. // #ifndef ThePEG_RandomGenerator_H #define ThePEG_RandomGenerator_H // This is the declaration of the RandomGenerator class. #include "ThePEG/Config/ThePEG.h" #include "ThePEG/Interface/Interfaced.h" #include "gsl/gsl_rng.h" namespace ThePEG { /** * RandomGenerator is an interface to the CLHEP::RandomEngine * classes. To avoid excessive virtual function calls, the * RandomGenerator caches random numbers generated by the engine which * are then retrieved by the non-virtual inlined rnd() method. * * Sub-classes of RandomGenerator should be used to * implement a particular random engine. * * RandomGenerator only provides a flat distribution between 0 and * 1. Any other distribution can be achieved using the CLHEP random * classes using the engine returned from the randomGenerator() * method. * * @see \ref RandomGeneratorInterfaces "The interfaces" * defined for RandomGenerator. * @see UseRandom */ class RandomGenerator: public Interfaced { public: /** A vector of doubles. */ typedef vector RndVector; /** The size_type of RndVector. */ typedef RndVector::size_type size_type; public: /** @name Standard constructors and destructors. */ //@{ /** * Default constructor. */ RandomGenerator(); /** * Copy-constructor. */ RandomGenerator(const RandomGenerator &); /** * Destructor. */ virtual ~RandomGenerator(); //@} /** * Reset the underlying random engine with the given \a seed. If the * \a seed is set to -1 a standard seed will be used. */ virtual void setSeed(long seed) = 0; /** @name Functions to return random numbers. */ //@{ /** * Return a (possibly cached) flat random number in the interval * \f$]0,1[\f$. */ double rnd() { if ( nextNumber == theNumbers.end() ) fill(); return *nextNumber++; } /** * Return a flat random number in the interval * \f$]0,b[\f$. */ template Unit rnd(Unit b) { return b*rnd(); } /** * Return a flat random number in the interval * \f$]a,b[\f$. */ template Unit rnd(Unit a, Unit b) { return a + rnd(b - a); } /** * Return \a n flat random number in the interval * \f$]0,1[\f$. */ RndVector rndvec(int n) { RndVector ret(n); for ( int i = 0; i < n; ++i ) ret[i] = rnd(); return ret; } /** * Return a (possibly cached) flat random number in the interval * \f$]0,1[\f$. */ double operator()() { return rnd(); } /** * Return a (possibly cached) flat integer random number in the * interval \f$[0,N[\f$. * Function was introduced since otherwise operator()() is used if a double is given * resulting in a \f$]0,1[\f$ distribution. */ double operator()(double N) { return double(rnd() * N); } /** * Return a (possibly cached) flat integer random number in the * interval \f$[0,N[\f$. */ long operator()(long N) { return long(rnd() * N); } /** * Return a true with probability \a p. Uses rnd(). */ bool rndbool(double p = 0.5) { return rnd() < p; } /** * Return a true with probability \a p. Uses rnd(). Also uses * push_back(double). */ bool prndbool(double p = 0.5); /** * Return a true with probability \a p1/(\a p1+\a p2). Uses * rnd(). */ bool rndbool(double p1, double p2) { return rndbool(p1/(p1 + p2)); } /** * Return a true with probability \a p1/(\a p1+\a p2). Uses * rnd(). Also uses push_back(double). */ bool prndbool(double p1, double p2) { return prndbool(p1/(p1 + p2)); } /** * Return -1, 0, or 1 with relative probabilities \a p1, \a p2, \a * p3. Uses rnd(). */ int rndsign(double p1, double p2, double p3); /** * Return -1, 0, or 1 with relative probabilities \a p1, \a p2, \a * p3. Uses rnd(). Also uses push_back(double). */ int prndsign(double p1, double p2, double p3); /** * Return an integer \f$i\f$ with probability p\f$i\f$/(\a p0+\a * p1). Uses rnd(). */ int rnd2(double p0, double p1) { return rndbool(p0, p1)? 0: 1; } /** * Return an integer \f$i\f$ with probability p\f$i\f$/(\a p0+\a * p1+\a p2). Uses rnd(). */ int rnd3(double p0, double p1, double p2) { return 1 + rndsign(p0, p1, p2); } /** * Return an integer/ \f$i\f$ with probability p\f$i\f$(\a p0+\a * p1+\a p2+\a p3). Uses rnd(). */ int rnd4(double p0, double p1, double p2, double p3); /** + * Return an integer/ \f$i\f$ with probability p\f$i\f$(\a p0+\a + * p1+\a p2+\a p3+\a p4). Uses rnd(). + */ + int rnd5(double p0, double p1, double p2, double p3, double p4); + + /** * Return a number between zero and infinity, distributed according * to \f$e^-x\f$. */ double rndExp() { return -log(rnd()); } /** * Return a number between zero and infinity, distributed according * to \f$e^-{x/\mu}\f$ where \f$\mu\f$ is the \a mean value. */ template Unit rndExp(Unit mean) { return mean*rndExp(); } /** * Return two numbers distributed according to a Gaussian distribution * with zero mean and unit variance. * * @param[out] First random number * @param[out] Second random number */ void rndGaussTwoNumbers(double & randomNumberOne, double & randomNumberTwo) { double r = sqrt(-2.0*log(rnd())); double phi = rnd()*2.0*Constants::pi; randomNumberOne = r*sin(phi); randomNumberTwo = r*cos(phi); } /** * Return a number distributed according to a Gaussian distribution * with zero mean and unit variance. * A second number is cached and returned the next time. * This function calls the rndGaussTwoNumbers function which returns two numbers at once. */ double rndGauss() { if ( gaussSaved ) { gaussSaved = false; return savedGauss; } double randomNumberOne, randomNumberTwo; rndGaussTwoNumbers(randomNumberOne, randomNumberTwo); savedGauss = randomNumberTwo; gaussSaved = true; return randomNumberOne; } /** * Return a number distributed according to a Gaussian distribution * with a given standard deviation, \a sigma, and a given \a mean. */ template Unit rndGauss(Unit sigma, Unit mean = Unit()) { return mean + sigma*rndGauss(); } /** * Return two numbers distributed according to a Gaussian distribution * with a given standard deviation, \a sigma, and a given \a mean. */ template void rndGaussTwoNumbers(Unit & randomNumberOne, Unit & randomNumberTwo, Unit sigma, Unit mean = Unit()) { double r1,r2; rndGaussTwoNumbers(r1,r2); randomNumberOne = mean + sigma * r1; randomNumberTwo = mean + sigma * r2; } /** * Return a positive number distributed according to a * non-relativistic Breit-Wigner with a given width, \a gamma, and a * given \a mean. */ template Unit rndBW(Unit mean, Unit gamma) { if ( gamma <= Unit() ) return mean; return mean + 0.5*gamma*tan(rnd(atan(-2.0*mean/gamma), Constants::pi/2)); } /** * Return a positive number distributed according to a * non-relativistic Breit-Wigner with a given width, \a gamma, and a * given \a mean. The distribution is cut-off so that the number is * between \a mean - \a cut and \a mean + \a cut */ template Unit rndBW(Unit mean, Unit gamma, Unit cut) { if ( gamma <= Unit() || cut <= Unit() ) return mean; return mean + 0.5*gamma*tan(rnd(atan(-2.0*min(mean,cut)/gamma), atan(2.0*cut/gamma))); } /** * Return a positive number distributed according to a relativistic * Breit-Wigner with a given width, \a gamma, and a given \a mean. */ template Unit rndRelBW(Unit mean, Unit gamma) { if ( gamma <= Unit() ) return mean; return sqrt(sqr(mean) + mean*gamma*tan(rnd(atan(-mean/gamma), Constants::pi/2))); } /** * Return a positive number distributed according to a relativistic * Breit-Wigner with a given width, \a gamma, and a given \a * mean. The distribution is cut-off so that the number is between * \a mean - \a cut and \a mean + \a cut */ template Unit rndRelBW(Unit mean, Unit gamma, Unit cut) { if ( gamma <= Unit() || cut <= Unit() ) return mean; double minarg = cut > mean? -mean/gamma: (sqr(mean - cut) - sqr(mean))/(gamma*mean); double maxarg = (sqr(mean + cut) - sqr(mean))/(mean*gamma); return sqrt(sqr(mean) + mean*gamma*tan(rnd(atan(minarg), atan(maxarg)))); } /** * Return a non-negative number generated according to a Poissonian * distribution with a given \a mean. Warning: the method * implemented is very slow for large mean and large return * values. For this reason the maximum return value is given by \a * nmax. */ long rndPoisson(double mean); //@} /** @name Access the cached random numbers from the underlying engine. */ //@{ /** * Give back a partly unused random number. This is typically used * when generating integral random numbers. In eg. rndbool(double p) * a random number r is drawn and if it is less than * p true is returned, but r is still a * good random number in the interval ]0,p[. Hence * r/p is still a good random number in the interval * ]0,1[ and this is then pushed back into the cache * and is used by the next call to rnd(). Note that the resulting * random number is of lesser quality, and successive calls to * push_back() should be avoided. To ensure a highest quality random * number random number in the next call to rnd(), pop_back() should * be used. */ void push_back(double r) { if ( r > 0.0 && r < 1.0 && nextNumber != theNumbers.begin() ) *--nextNumber = r; } /** * Discard the next random number in the cache. */ void pop_back() { if ( nextNumber != theNumbers.end() ) ++nextNumber; } /** * Discard all random numbers in the cache. Typically used after the * internal random engine has been reinitialized for some reason. */ void flush() { nextNumber = theNumbers.end(); gaussSaved = false; } /** * Generate n random numbers between 0 and 1 and put them in the * output iterator. */ template void rnd(OutputIterator o, size_type n) { while ( n-- ) *o++ = rnd(); } //@} protected: /** * Initializes this random generator. This should be done first of * all before the initialization of any other object associated with * an event generator. */ virtual void doinit(); public: /** @name Functions used by the persistent I/O system. */ //@{ /** * Function used to write out object persistently. * @param os the persistent output stream written to. */ void persistentOutput(PersistentOStream & os) const; /** * Function used to read in object persistently. * @param is the persistent input stream read from. * @param version the version number of the object when written. */ void persistentInput(PersistentIStream & is, int version); //@} /** * Standard Init function used to initialize the interface. */ static void Init(); /** * Return a gsl_rng interface to this random generator. */ gsl_rng * getGslInterface() { return gsl; } protected: /** * Utility function for the interface. */ void setSize(size_type newSize); /** * Fill the cache with random numbers. */ virtual void fill() = 0; /** * A vector of cached numbers. */ RndVector theNumbers; /** * Iterator pointing to the next number to be extracted */ RndVector::iterator nextNumber; /** * The size of the cache. */ size_type theSize; /** * The seed to initialize the random generator with. */ long theSeed; /** * A saved Gaussian random number. */ mutable double savedGauss; /** * Indicate the precense of a saved Gaussian random number. */ mutable bool gaussSaved; /** * A pinter to a gsl_rng interface to this generator. */ gsl_rng * gsl; private: /** * Describe a concrete class with persistent data. Note that the * class should in principle be abstract. */ static ClassDescription initRandomGenerator; /** * Private and non-existent assignment operator. */ RandomGenerator & operator=(const RandomGenerator &) = delete; }; /** @cond TRAITSPECIALIZATIONS */ /** This template specialization informs ThePEG about the base classes * of RandomGenerator. */ template <> struct BaseClassTrait: public ClassTraitsType { /** Typedef of the first base class of RandomGenerator. */ typedef Interfaced NthBase; }; /** This template specialization informs ThePEG about the name of the * RandomGenerator class. */ template <> struct ClassTraits: public ClassTraitsBase { /** Return a platform-independent class name */ static string className() { return "ThePEG::RandomGenerator"; } /** This class should in principle be abstract, therefore the create() method will throw a std::logic_error if called. */ static TPtr create() { throw std::logic_error("Tried to instantiate abstract class" "'Pythis7::RandomGenerator'"); } }; /** @endcond */ } #endif /* ThePEG_RandomGenerator_H */ diff --git a/Repository/UseRandom.h b/Repository/UseRandom.h --- a/Repository/UseRandom.h +++ b/Repository/UseRandom.h @@ -1,289 +1,297 @@ // -*- C++ -*- // // UseRandom.h is a part of ThePEG - Toolkit for HEP Event Generation // Copyright (C) 1999-2017 Leif Lonnblad // // ThePEG is licenced under version 3 of the GPL, see COPYING for details. // Please respect the MCnet academic guidelines, see GUIDELINES for details. // #ifndef ThePEG_UseRandom_H #define ThePEG_UseRandom_H // This is the declaration of the UseRandom class. #include "ThePEG/Repository/RandomGenerator.h" namespace ThePEG { /** * This UseRandom class keeps a static stack of RandomGenerator * objects which can be used anywhere by any class. When an * EventGenerator is initialized or run it adds a RandomGenerator * object to the stack which can be used by any other object being * initialized or run through the static functions of the UseRandom * class. If someone needs to use an alternative RandomGenerator * object a new UseRandom object can be constructed with a pointer to * the desired RandomGenerator object as argument and that object will * the be used by the static UseRandom functions until the UseRandom * object is destructed. * * @see RandomGenerator * @see EventGenerator * */ class UseRandom { public: /** * Default constructor does nothing. */ UseRandom() : randomPushed(false) {} /** * Copy-constructor does nothing. */ UseRandom(const UseRandom &) : randomPushed(false) {} /** * Construct a new object specifying a new RandomGenerator, \a r, to * be used during this objects lifetime */ UseRandom(const RanGenPtr & r) : randomPushed(false) { if ( r ) { theRandomStack.push_back(r); randomPushed = true; } } /** * The destructor removing the RandomGenerator specified in the * constructor from the stack. */ ~UseRandom() { if ( randomPushed ) theRandomStack.pop_back(); } public: /** * Return a reference to the currently chosen RandomGenerator object. */ static RandomGenerator & current() { return *theRandomStack.back(); } /** * Return a pointer to the currently chosen RandomGenerator object. */ // static RandomEngine * currentEngine() { // return &(current().randomGenerator()); // } /** * Return a simple flat random number (from the current * RandomGenerator object) in the range ]0,1[. */ static double rnd() { return current().rnd(); } /** * Return \a n simple flat random number (from the current * RandomGenerator object) in the range ]0,1[. */ static RandomGenerator::RndVector rndvec(int n) { return current().rndvec(n); } /** * Return a simple flat random number (from the current * RandomGenerator object) in the range ]0,\a xu[. */ template static Unit rnd(Unit xu) { return current().rnd(xu); } /** * Return a simple flat random number (from the current * RandomGenerator object) in the range ]\a xl,\a xu[. */ template static Unit rnd(Unit xl, Unit xu) { return current().rnd(xl, xu); } /** * Return a true with probability \a p (default 0.5). */ static bool rndbool(double p = 0.5) { return current().rndbool(p); } /** * Return a true with probability \a p (default 0.5). Uses push_back * to reuse random number. */ static bool prndbool(double p = 0.5) { return current().rndbool(p); } /** * Return a true with probability \a p1/(\a p1+\a p2). */ static bool rndbool(double p1, double p2) { return current().rndbool(p1, p2); } /** * Return a true with probability \a p1/(\a p1+\a p2). Uses * push_back to reuse random number. */ static bool prndbool(double p1, double p2) { return current().rndbool(p1, p2); } /** * Return -1, 0, or 1 with relative probabilities \a p1, \a p2, \a * p3. */ static int rndsign(double p1, double p2, double p3) { return current().rndsign(p1, p2, p3); } /** * Return -1, 0, or 1 with relative probabilities \a p1, \a p2, \a * p3. Uses push_back to reuse random number. */ static int prndsign(double p1, double p2, double p3) { return current().rndsign(p1, p2, p3); } /** * Return an integer \f$i\f$ with probability p\f$i\f$/(\a p0+\a * p1). */ static int rnd2(double p0, double p1) { return current().rnd2(p0, p1); } /** * Return an integer \f$i\f$ with probability p\f$i\f$/(\a p0+\a * p1+\a p2). */ static int rnd3(double p0, double p1, double p2) { return current().rnd3(p0, p1, p2); } /** * Return an integer/ \f$i\f$ with probability p\f$i\f$(\a p0+\a * p1+\a p2+\a p3). */ static int rnd4(double p0, double p1, double p2, double p3) { return current().rnd4(p0, p1, p2, p3); } /** + * Return an integer/ \f$i\f$ with probability p\f$i\f$(\a p0+\a + * p1+\a p2+\a p3+\a p4). + */ + static int rnd5(double p0, double p1, double p2, double p3, double p4) { + return current().rnd5(p0, p1, p2, p3, p4); + } + + /** * Return a simple flat random integrer number in the range [0,\a xu[. */ static long irnd(long xu = 2) { return long(rnd() * xu); } /** * Return a simple flat random integrer number in the range [\a xl,\a xu[. */ static long irnd(long xl, long xu) { return xl + irnd(xu-xl); } /** * Return a number between zero and infinity, distributed according * to \f$e^-x\f$. */ static double rndExp() { return current().rndExp(); } /** * Return a number between zero and infinity, distributed according * to \f$e^-{x/\mu}\f$ where \f$\mu\f$ is the \a mean value. */ template static Unit rndExp(Unit mean) { return current().rndExp(mean); } /** * Return a number distributed according to a Gaussian distribution * with zero mean and unit variance. */ static double rndGauss() { return current().rndGauss(); } /** * Return a number distributed according to a Gaussian distribution * with a given standard deviation, \a sigma, and a given \a mean. */ template static Unit rndGauss(Unit sigma, Unit mean = Unit()) { return current().rndGauss(sigma, mean); } /** * Return a positive number distributed according to a * non-relativistic Breit-Wigner with a given width, \a gamma, and a * given \a mean. */ template static Unit rndBW(Unit mean, Unit gamma) { return current().rndBW(mean, gamma); } /** * Return a positive number distributed according to a * non-relativistic Breit-Wigner with a given width, \a gamma, and a * given \a mean. The distribution is cut-off so that the number is * between \a mean - \a cut and \a mean + \a cut */ template static Unit rndBW(Unit mean, Unit gamma, Unit cut) { return current().rndBW(mean, gamma, cut); } /** * Return a positive number distributed according to a relativistic * Breit-Wigner with a given width, \a gamma, and a given \a mean. */ template static Unit rndRelBW(Unit mean, Unit gamma) { return current().rndRelBW(mean, gamma); } /** * Return a positive number distributed according to a relativistic * Breit-Wigner with a given width, \a gamma, and a given \a * mean. The distribution is cut-off so that the number is between * \a mean - \a cut and \a mean + \a cut */ template static Unit rndRelBW(Unit mean, Unit gamma, Unit cut) { return current().rndRelBW(mean, gamma, cut); } /** * Return a non-negative number generated according to a Poissonian * distribution with a given \a mean. */ static long rndPoisson(double mean) { return current().rndPoisson(mean); } private: /** * The stack of RandomGenerators requested. */ static vector theRandomStack; /** * True if this object is responsible for pushing a RandomGenerator * onto the stack. */ bool randomPushed; private: /** * Private and non-existent assignment operator. */ UseRandom & operator=(const UseRandom &) = delete; }; } #endif /* ThePEG_UseRandom_H */