Ignore:
Timestamp:
Dec 7, 2005, 11:01:17 AM (16 years ago)
Author:
Jari Häkkinen
Message:

Fixed documentation and moved around code.

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/lib/random/random.h

    r411 r424  
    1212
    1313namespace theplu {
    14 namespace random { 
     14namespace random {
    1515
    1616  ///
    1717  /// @brief Random Number Generator
    1818  ///
    19   /// The RNG class is wrapper to the GSL random number generator
    20   /// (rng). This class provides a single global instance of the rng,
    21   /// and makes sure there is only one point of access to the
    22   /// generator.
     19  /// The RNG class is wrapper to the GSL random number generator
     20  /// (rng). This class provides a single global instance of the rng,
     21  /// and makes sure there is only one point of access to the
     22  /// generator.
    2323  ///
    2424  /// There is information about how to change seeding and generators
     
    4444  /// thread safety.
    4545  ///
    46   /// @see Design Patterns (the singleton and adapter pattern). GSL
    47   /// documentation.
     46  /// @see Design Patterns (the singleton and adapter pattern). GSL
     47  /// documentation.
     48  ///
     49  /// @todo Is this class properly implemented? The underlying random
     50  /// number genereator should be a singleton, while allowing for
     51  /// several different distribution to be used. Jari's feeling is
     52  /// that the current implementation restricts the use to only one
     53  /// distribution per binary!
    4854  ///
    4955  class RNG
     
    5359    virtual ~RNG(void);
    5460
    55     static RNG* instance(u_long seed=0);
     61    ///
     62    /// @brief Get an instance of the random number
     63    /// generator/distribution.
     64    ///
     65    /// @return A pointer to the random number generator.
     66    ///
     67    static RNG* instance(u_long seed=0)
     68      { if (!instance_) instance_=new RNG(seed); return instance_; }
    5669
    5770    ///
     
    8194    /// @brief Set the seed \a s for the rng.
    8295    ///
    83     /// @see seed_dev_urandom
     96    /// @see seed_from_devurandom
    8497    ///
    8598    inline void seed(u_long s) const { gsl_rng_set(rng_,s); }
     
    93106
    94107  private:
    95 
    96108    RNG(u_long seed);
    97109
     
    100112  };
    101113
    102 
    103 
    104   ///
    105   /// @brief Continuous random number distributions.
    106   ///
    107   /// Abstract base class for continuous random number distributions.
    108   ///
    109   class Continuous
    110   {
    111   public:
    112 
    113     ///
    114     /// @brief Constructor
    115     ///
    116     inline Continuous(void) { rng_=RNG::instance(); }
    117 
    118     ///
    119     /// @return A random number
    120     ///
    121     virtual double operator()(void) const = 0;
    122 
     114  // --------------------- Discrete distribtuions ---------------------
     115
     116  ///
     117  /// @brief Discrete random number distributions.
     118  ///
     119  /// Abstract base class for discrete random number
     120  /// distributions. Given K discrete events with different
     121  /// probabilities \f$ P[k] \f$, produce a random value k consistent
     122  /// with its probability.
     123  ///
     124  class Discrete
     125  {
     126  public:
     127    ///
     128    /// @brief Constructor
     129    ///
     130    inline Discrete(void) { rng_=RNG::instance(); }
     131
     132    ///
     133    /// @return A random number.
     134    ///
     135    virtual u_long operator()(void) const = 0;
     136   
    123137  protected:
    124138    RNG* rng_;
    125139  };
    126140
     141  ///
     142  /// @brief General
     143  ///
     144  class DiscreteGeneral : public Discrete
     145  {
     146  public:
     147    ///
     148    /// @brief Constructor
     149    ///
     150    /// @param hist is a Histogram defining the probability distribution
     151    ///
     152    DiscreteGeneral(const statistics::Histogram& hist);
     153   
     154    ///
     155    /// @brief Destructor
     156    ///
     157    virtual ~DiscreteGeneral(void);
     158
     159    ///
     160    /// The generated number is an integer and proportinal to the
     161    /// frequency in the corresponding histogram bin. In other words,
     162    /// the probability that 0 is returned is proportinal to the size
     163    /// of the first bin.
     164    ///
     165    /// @return A random number.
     166    ///
     167    inline u_long
     168    operator()(void) const { return gsl_ran_discrete(rng_->rng(), gen_); }
     169
     170  private:
     171     gsl_ran_discrete_t* gen_;
     172  };
     173
     174  ///
     175  /// @brief Discrete uniform distribution
     176  ///
     177  /// Discrete uniform distribution also known as the "equally likely
     178  /// outcomes" distribution. Each outcome, in this case an integer
     179  /// from [0,n-1] , have equal probability to occur.
     180  ///
     181  /// Distribution function \f$ p(k) = \frac{1}{n+1} \f$ for \f$ 0 \le
     182  /// k < n \f$ \n
     183  /// Expectation value: \f$ \frac{n-1}{2} \f$ \n
     184  /// Variance: \f$ \frac{1}{12}(n-1)(n+1) \f$
     185  ///
     186  class DiscreteUniform : public Discrete
     187  {
     188  public:
     189    ///
     190    /// @brief Default constructor.
     191    ///
     192    DiscreteUniform(void) : range_(rng_->max()+1) {}
     193
     194    ///
     195    /// @brief Constructor.
     196    ///
     197    /// @param n sets the range. Object will generate integers from
     198    /// [0,n-1].
     199    ///
     200    DiscreteUniform(const u_long n) : range_(n)
     201    { if ( range_-1>rng_->max() ) range_=rng_->max()+1; }
     202
     203    ///
     204    /// This function returns a random integer from 0 to n-1
     205    /// inclusive. All integers in the range [0,n-1] are equally
     206    /// likely. n is set in constructor.
     207    ///
     208    inline u_long
     209    operator()(void) const { return gsl_rng_uniform_int(rng_->rng(), range_); }
     210
     211    ///
     212    /// This function returns a random integer from 0 to n-1
     213    /// inclusive. All integers in the range [0,n-1] are equally
     214    /// likely.
     215    ///
     216    inline u_long operator()(const u_long n) const
     217    { return gsl_rng_uniform_int(rng_->rng(), n); }
     218
     219  private:
     220    u_long range_;
     221  };
     222
     223  ///
     224  /// @brief Poisson Distribution
     225  ///
     226  /// Having a Poisson process (i.e. no memory), number of occurences
     227  /// within a given time window is Poisson distributed. This
     228  /// distribution is the limit of a Binomial distribution when number
     229  /// of attempts is large, and the probability for one attempt to be
     230  /// succesful is small (in such a way that the expected number of
     231  /// succesful attempts is \f$ m \f$.
     232  ///
     233  /// Probability function \f$ p(k) = e^{-m}\frac{m^k}{k!} \f$ for \f$ 0 \le
     234  /// k  \f$ \n
     235  /// Expectation value: \f$ m \f$ \n
     236  /// Variance: \f$ m \f$
     237  ///
     238  class Poisson : public Discrete
     239  {
     240  public:
     241    ///
     242    /// @brief Constructor
     243    ///
     244    /// @param m is expectation value
     245    inline Poisson(const double m=1) : m_(m) {}
     246
     247    ///
     248    /// @return A Poisson distributed number.
     249    ///
     250    inline u_long
     251    operator()(void) const { return gsl_ran_poisson(rng_->rng(), m_); }
     252
     253    ///
     254    /// @return A Poisson distributed number with expectation value \a
     255    /// m
     256    ///
     257    /// @note this operator ignores parameters set in Constructor
     258    ///
     259    inline u_long
     260    operator()(const double m) const { return gsl_ran_poisson(rng_->rng(), m); }
     261
     262  private:
     263    double m_;
     264  };
     265
     266  // --------------------- Continuous distribtuions ---------------------
     267
     268  ///
     269  /// @brief Continuous random number distributions.
     270  ///
     271  /// Abstract base class for continuous random number distributions.
     272  ///
     273  class Continuous
     274  {
     275  public:
     276
     277    ///
     278    /// @brief Constructor
     279    ///
     280    inline Continuous(void) { rng_=RNG::instance(); }
     281
     282    ///
     283    /// @return A random number
     284    ///
     285    virtual double operator()(void) const = 0;
     286
     287  protected:
     288    RNG* rng_;
     289  };
    127290
    128291  ///
     
    139302  {
    140303  public:
    141 
    142304    inline double operator()(void) const { return gsl_rng_uniform(rng_->rng());}
    143 
     305  };
     306
     307  ///
     308  /// Class to generate numbers from a histogram in a continuous manner.
     309  ///
     310  class ContinuousGeneral : public Continuous
     311  {
     312  public:
     313    ///
     314    /// @brief Constructor
     315    ///
     316    /// @param hist is a Histogram defining the probability distribution
     317    ///
     318    inline ContinuousGeneral(const statistics::Histogram& hist)
     319      : discrete_(DiscreteGeneral(hist)), hist_(hist) {}
     320
     321    ///
     322    /// @brief Destructor
     323    ///
     324    virtual ~ContinuousGeneral(void);
     325
     326    ///
     327    /// The number is generated in a two step process. First the bin
     328    /// in the histogram is randomly selected (see
     329    /// DiscreteGeneral). Then a number is generated uniformly from
     330    /// the interval defined by the bin.
     331    ///
     332    /// @return A random number.
     333    ///
     334    inline double operator()(void) const
     335    { return hist_.observation_value(discrete_())+(u_()-0.5)*hist_.spacing(); }
     336
     337  private:
     338    const DiscreteGeneral discrete_;
     339    const statistics::Histogram& hist_;
     340    ContinuousUniform u_;
     341  };
     342
     343  ///
     344  /// @brief Generator of random numbers from an exponential
     345  /// distribution.
     346  ///
     347  /// The distribution function is \f$ f(x) = \frac{1}{m}\exp(-x/a)
     348  /// \f$ for \f$ x \f$ with the expectation value \f$ m \f$ and
     349  /// variance \f$ m^2 \f$
     350  ///
     351  class Exponential : public Continuous
     352  {
     353  public:
     354    ///
     355    /// @brief Constructor
     356    ///
     357    /// @param \a m is the expectation value of the distribution.
     358    ///
     359    inline Exponential(const double m=1) : m_(m) {}
     360
     361    ///
     362    /// @return A random number from exponential distribution.
     363    ///
     364    inline double
     365    operator()(void) const { return gsl_ran_exponential(rng_->rng(), m_); }
     366
     367    ///
     368    /// @return A random number from exponential distribution, with
     369    /// expectation value \a m
     370    ///
     371    /// @note This operator ignores parameters given in constructor.
     372    ///
     373    inline double operator()(const double m) const
     374    { return gsl_ran_exponential(rng_->rng(), m); }
     375
     376  private:
     377    double m_;
    144378  };
    145379
     
    165399    /// m is the expectation value \f$ \mu \f$ of the distribution
    166400    ///
    167     inline Gaussian(const double s=1, const double m=0)
    168       : m_(m), s_(s) {}
     401    inline Gaussian(const double s=1, const double m=0) : m_(m), s_(s) {}
    169402
    170403    ///
    171404    /// @return A random Gaussian number
    172405    ///
    173     inline double operator()(void) const
    174     { return gsl_ran_gaussian(rng_->rng(), s_)+m_; }
     406    inline double
     407    operator()(void) const { return gsl_ran_gaussian(rng_->rng(), s_)+m_; }
    175408
    176409    ///
     
    179412    /// @note this operator ignores parameters given in Constructor
    180413    ///
    181     inline double operator()(const double s) const
    182     { return gsl_ran_gaussian(rng_->rng(), s); }
     414    inline double
     415    operator()(const double s) const { return gsl_ran_gaussian(rng_->rng(), s); }
    183416
    184417    ///
     
    187420    /// @note this operator ignores parameters given in Constructor
    188421    ///
    189     inline double operator()(const double s, const double m) const 
     422    inline double operator()(const double s, const double m) const
    190423    { return gsl_ran_gaussian(rng_->rng(), s)+m; }
    191424
     
    195428  };
    196429
    197   ///
    198   /// @brief Exponential distribution
    199   ///
    200   /// Class for generating a random number from an Exponential
    201   /// distribution.
    202   ///
    203   /// Distribution function \f$ f(x) = \frac{1}{m}\exp(-x/a) \f$ for
    204   /// \f$ x \f$ \n
    205   /// Expectation value: \f$ m \f$ \n
    206   /// Variance: \f$ m^2 \f$
    207   ///
    208   class Exponential : public Continuous
    209   {
    210   public:
    211     ///
    212     /// @brief Constructor
    213     /// @param m is the expectation value of the distribution.
    214     ///
    215     inline Exponential(const double m=1)
    216       : m_(m){}
    217 
    218     ///
    219     /// @return A random number from exponential distribution
    220     ///
    221     inline double operator()(void) const
    222     { return gsl_ran_exponential(rng_->rng(), m_); }
    223 
    224     ///
    225     /// @return A random number from exponential distribution, with
    226     /// expectation value \a m
    227     /// @note this operator ignores parameters given in Constructor
    228     ///
    229     inline double operator()(const double m) const
    230     { return gsl_ran_exponential(rng_->rng(), m); }
    231 
    232   private:
    233     double m_;
    234   };
    235 
    236   ///
    237   /// @brief Discrete random number distributions.
    238   ///
    239   /// Abstract base class for discrete random number
    240   /// distributions. Given K discrete events with different
    241   /// probabilities \f$ P[k] \f$, produce a random value k consistent
    242   /// with its probability.
    243   ///
    244   class Discrete
    245   {
    246   public:
    247     ///
    248     /// @brief Constructor
    249     ///
    250     inline Discrete(void) { rng_=RNG::instance(); }
    251 
    252     ///
    253     /// @return A random number.
    254     ///
    255     virtual u_long operator()(void) const = 0;
    256    
    257   protected:
    258     RNG* rng_;
    259   };
    260 
    261   ///
    262   /// @brief Discrete uniform distribution
    263   ///
    264   /// Discrete uniform distribution also known as the "equally likely
    265   /// outcomes" distribution. Each outcome, in this case an integer
    266   /// from [0,n-1] , have equal probability to occur.
    267   ///
    268   /// Distribution function \f$ p(k) = \frac{1}{n+1} \f$ for \f$ 0 \le
    269   /// k < n \f$ \n
    270   /// Expectation value: \f$ \frac{n-1}{2} \f$ \n
    271   /// Variance: \f$ \frac{1}{12}(n-1)(n+1) \f$
    272   ///
    273 
    274   class DiscreteUniform : public Discrete
    275   {
    276   public:
    277     ///
    278     /// @brief Constructor.
    279     ///
    280     DiscreteUniform(void);
    281 
    282     ///
    283     /// @brief Constructor.
    284     ///
    285     /// @param n sets the range. Object will generate integers from
    286     /// [0,n-1].
    287     ///
    288     DiscreteUniform(const u_long n);
    289      
    290     ///
    291     /// This function returns a random integer from 0 to n-1
    292     /// inclusive. All integers in the range [0,n-1] are equally
    293     /// likely. n is set in constructor.
    294     ///
    295     inline u_long operator()(void) const
    296     { return gsl_rng_uniform_int(rng_->rng(), range_); }
    297 
    298     ///
    299     /// This function returns a random integer from 0 to n-1
    300     /// inclusive. All integers in the range [0,n-1] are equally
    301     /// likely.
    302     ///
    303     inline u_long operator()(const u_long n) const
    304     { return gsl_rng_uniform_int(rng_->rng(), n); }
    305 
    306 
    307   private:
    308     u_long range_;
    309   };
    310 
    311 
    312   ///
    313   /// @brief Poisson Distribution
    314   ///
    315   /// Having a Poisson process (no memory), number of occurences
    316   /// within a given time window is Poisson distributed. This
    317   /// distribution is the limit of a Binomial distribution when number
    318   /// of attempts is large, and the probability for one attempt to be
    319   /// succesful is small (in such a way that the expected number of
    320   /// succesful attempts is \f$ m \f$.
    321   ///
    322   /// Probability function \f$ p(k) = e^{-m}\frac{m^k}{k!} \f$ for \f$ 0 \le
    323   /// k  \f$ \n
    324   /// Expectation value: \f$ m \f$ \n
    325   /// Variance: \f$ m \f$
    326   ///
    327   class Poisson : public Discrete
    328   {
    329   public:
    330     ///
    331     /// @brief Constructor
    332     ///
    333     /// @param m is expectation value
    334     inline Poisson(const double m=1)
    335       : m_(m){}
    336 
    337     ///
    338     /// @return A Poisson distributed number.
    339     ///
    340     inline u_long operator()(void) const
    341     { return gsl_ran_poisson(rng_->rng(), m_); }
    342 
    343     ///
    344     /// @return A Poisson distributed number with expectation value \a
    345     /// m
    346     /// @note this operator ignores parameters set in Constructor
    347     inline u_long operator()(const double m) const
    348     { return gsl_ran_poisson(rng_->rng(), m); }
    349 
    350   private:
    351     double m_;
    352   };
    353 
    354   ///
    355   /// @brief General
    356   ///
    357   class DiscreteGeneral : public Discrete
    358   {
    359   public:
    360     ///
    361     /// @brief Constructor
    362     ///
    363     /// @param hist is a Histogram defining the probability distribution
    364     ///
    365     DiscreteGeneral(const statistics::Histogram& hist) ;
    366    
    367     ///
    368     /// @brief Destructor
    369     ///
    370     virtual ~DiscreteGeneral(void);
    371 
    372     ///
    373     /// The generated number is an integer and proportinal to the
    374     /// frequency in the corresponding histogram bin. In other words,
    375     /// the probability that 0 is returned is proportinal to the size
    376     /// of the first bin.
    377     ///
    378     /// @return A random number.
    379     ///
    380     inline u_long operator()(void) const
    381     { return gsl_ran_discrete(rng_->rng(), gen_); }
    382 
    383   private:
    384      gsl_ran_discrete_t* gen_;
    385   };
    386 
    387 
    388   ///
    389   /// Class to generate numbers from a histogram in a continuous manner.
    390   ///
    391   class ContinuousGeneral : public Continuous
    392   {
    393   public:
    394     ///
    395     /// @brief Constructor
    396     ///
    397     /// @param hist is a Histogram defining the probability distribution
    398     ///
    399     inline ContinuousGeneral(const statistics::Histogram& hist)
    400       : discrete_(DiscreteGeneral(hist)), hist_(hist) {}
    401    
    402     ///
    403     /// @brief Destructor
    404     ///
    405     virtual ~ContinuousGeneral(void);
    406 
    407     ///
    408     /// The number is generated in a two step process. First the bin
    409     /// in the histogram is randomly selected (see
    410     /// DiscreteGeneral). Then a number is generated uniformly from
    411     /// the interval defined by the bin.
    412     ///
    413     /// @return A random number.
    414     ///
    415     inline double operator()(void) const
    416     { return hist_.observation_value(discrete_())+(u_()-0.5)*hist_.spacing(); }
    417 
    418   private:
    419     const DiscreteGeneral discrete_;
    420     const statistics::Histogram& hist_;
    421     ContinuousUniform u_;
    422   };
    423 
    424430
    425431}} // of namespace random and namespace theplu
Note: See TracChangeset for help on using the changeset viewer.