Statistics
| Branch: | Revision:

root / include / cnumgen.h @ f857128b

History | View | Annotate | Download (13.4 KB)

1
//==========================================================================
2
// DISTRIB.H
3
//
4
//                     OMNeT++/OMNEST
5
//            Discrete System Simulation in C++
6
//
7
// Random variate generation class
8
//
9
//
10
//==========================================================================
11

    
12

    
13
#ifndef __CNUMGEN_H_
14
#define __CNUMGEN_H_
15

    
16
#include "simkerneldefs.h"
17
#include "random.h"
18

    
19
class cMersenneTwister;
20

    
21
class SIM_API cNumberGenerator // : public cObject
22
{
23
  private:
24

    
25
        int nrRNGs;
26

    
27
        cRNG** rngs;
28

    
29
    /**
30
     * internal
31
     */
32
    double gamma_Marsaglia2000(double a, int rng);
33

    
34
    /**
35
     * internal
36
     */
37
    double gamma_MarsagliaTransf(double alpha, int rng);
38

    
39
    /**
40
     * internal
41
     */
42
    double unit_normal(int rng);
43

    
44
    /**
45
         * internal
46
         */
47
        void checkBounds(int rngId);
48

    
49
  public:
50

    
51
    /**
52
     * constructor
53
     */
54
    cNumberGenerator();
55

    
56
    /**
57
     * destructor
58
     */
59
    virtual ~cNumberGenerator();
60

    
61
    /*
62
     * Sets up seeds and instanciates the RNGs
63
     */
64
    void setupSeeds(unsigned nrRng, unsigned long * seeds);
65

    
66
    /**
67
     * @name Continuous distributions
68
     *
69
     * @ingroup RandomNumbers
70
     */
71
    //@{
72
    /**
73
     * Returns a random variate with uniform distribution in the range [a,b).
74
     *
75
     * @param a, b the interval, a<b
76
     * @param rng the underlying random number generator
77
     */
78
    double uniform(double a, double b, int =0);
79

    
80
    /**
81
     * Returns a random variate from the exponential distribution with the
82
     * given mean (that is, with parameter lambda=1/mean).
83
     *
84
     * @param mean mean value
85
     * @param rng the underlying random number generator
86
     */
87
    double exponential(double mean, int rng=0);
88

    
89
    /**
90
     * Returns a random variate from the normal distribution with the given mean
91
     * and standard deviation.
92
     *
93
     * @param mean mean of the normal distribution
94
     * @param stddev standard deviation of the normal distribution
95
     * @param rng the underlying random number generator
96
     */
97
    double normal(double mean, double stddev, int rng=0);
98

    
99
    /**
100
     * Normal distribution truncated to nonnegative values.
101
     * It is implemented with a loop that discards negative values until
102
     * a nonnegative one comes. This means that the execution time is not bounded:
103
     * a large negative mean with much smaller stddev is likely to result
104
     * in a large number of iterations.
105
     *
106
     * The mean and stddev parameters serve as parameters to the normal
107
     * distribution <i>before</i> truncation. The actual random variate returned
108
     * will have a different mean and standard deviation.
109
     *
110
     * @param mean mean of the normal distribution
111
     * @param stddev standard deviation of the normal distribution
112
     * @param rng the underlying random number generator
113
     */
114
    double truncnormal(double mean, double stddev, int rng=0);
115

    
116
    /**
117
     * Returns a random variate from the gamma distribution with parameters
118
     * alpha>0, theta>0. Alpha is known as the "shape" parameter, and theta
119
     * as the "scale" parameter.
120
     *
121
     * Some sources in the literature use the inverse scale parameter
122
     * beta = 1 / theta, called the "rate" parameter. Various other notations
123
     * can be found in the literature; our usage of (alpha,theta) is consistent
124
     * with Wikipedia and Mathematica (Wolfram Research).
125
     *
126
     * Gamma is the generalization of the Erlang distribution for non-integer
127
     * k values, which becomes the alpha parameter. The chi-square distribution
128
     * is a special case of the gamma distribution.
129
     *
130
     * For alpha=1, Gamma becomes the exponential distribution with mean=theta.
131
     *
132
     * The mean of this distribution is alpha*theta, and variance is alpha*theta<sup>2</sup>.
133
     *
134
     * Generation: if alpha=1, it is generated as exponential(theta).
135
     *
136
     * For alpha>1, we make use of the acceptance-rejection method in
137
     * "A Simple Method for Generating Gamma Variables", George Marsaglia and
138
     * Wai Wan Tsang, ACM Transactions on Mathematical Software, Vol. 26, No. 3,
139
     * September 2000.
140
     *
141
     * The alpha<1 case makes use of the alpha>1 algorithm, as suggested by the
142
     * above paper.
143
     *
144
     * @remark the name gamma_d is chosen to avoid ambiguity with
145
     * a function of the same name
146
     *
147
     * @param alpha >0  the "shape" parameter
148
     * @param theta >0  the "scale" parameter
149
     * @param rng the underlying random number generator
150
     */
151
    double gamma_d(double alpha, double theta, int rng=0);
152

    
153
    /**
154
     * Returns a random variate from the beta distribution with parameters
155
     * alpha1, alpha2.
156
     *
157
     * Generation is using relationship to Gamma distribution: if Y1 has gamma
158
     * distribution with alpha=alpha1 and beta=1 and Y2 has gamma distribution
159
     * with alpha=alpha2 and beta=2, then Y = Y1/(Y1+Y2) has beta distribution
160
     * with parameters alpha1 and alpha2.
161
     *
162
     * @param alpha1, alpha2 >0
163
     * @param rng the underlying random number generator
164
     */
165
    double beta(double alpha1, double alpha2, int rng=0);
166

    
167
    /**
168
     * Returns a random variate from the Erlang distribution with k phases
169
     * and mean mean.
170
     *
171
     * This is the sum of k mutually independent random variables, each with
172
     * exponential distribution. Thus, the kth arrival time
173
     * in the Poisson process follows the Erlang distribution.
174
     *
175
     * Erlang with parameters m and k is gamma-distributed with alpha=k
176
     * and beta=m/k.
177
     *
178
     * Generation makes use of the fact that exponential distributions
179
     * sum up to Erlang.
180
     *
181
     * @param k number of phases, k>0
182
     * @param mean >0
183
     * @param rng the underlying random number generator
184
     */
185
    double erlang_k(unsigned int k, double mean, int rng=0);
186

    
187
    /**
188
     * Returns a random variate from the chi-square distribution
189
     * with k degrees of freedom.  The chi-square distribution arises
190
     * in statistics. If Yi are k independent random variates from the normal
191
     * distribution with unit variance, then the sum-of-squares (sum(Yi^2))
192
     * has a chi-square distribution with k degrees of freedom.
193
     *
194
     * The expected value of this distribution is k. Chi_square with parameter
195
     * k is gamma-distributed with alpha=k/2, beta=2.
196
     *
197
     * Generation is using relationship to gamma distribution.
198
     *
199
     * @param k degrees of freedom, k>0
200
     * @param rng the underlying random number generator
201
     */
202
    double chi_square(unsigned int k, int rng=0);
203

    
204
    /**
205
     * Returns a random variate from the student-t distribution with
206
     * i degrees of freedom. If Y1 has a normal distribution and Y2 has a chi-square
207
     * distribution with k degrees of freedom then X = Y1 / sqrt(Y2/k)
208
     * has a student-t distribution with k degrees of freedom.
209
     *
210
     * Generation is using relationship to gamma and chi-square.
211
     *
212
     * @param i degrees of freedom, i>0
213
     * @param rng the underlying random number generator
214
     */
215
    double student_t(unsigned int i, int rng=0);
216

    
217
    /**
218
     * Returns a random variate from the Cauchy distribution (also called
219
     * Lorentzian distribution) with parameters a,b where b>0.
220
     *
221
     * This is a continuous distribution describing resonance behavior.
222
     * It also describes the distribution of horizontal distances at which
223
     * a line segment tilted at a random angle cuts the x-axis.
224
     *
225
     * Generation uses inverse transform.
226
     *
227
     * @param a
228
     * @param b  b>0
229
     * @param rng the underlying random number generator
230
     */
231
    double cauchy(double a, double b, int rng=0);
232

    
233
    /**
234
     * Returns a random variate from the triangular distribution with parameters
235
     * a <= b <= c.
236
     *
237
     * Generation uses inverse transform.
238
     *
239
     * @param a, b, c   a <= b <= c
240
     * @param rng the underlying random number generator
241
     */
242
    double triang(double a, double b, double c, int rng=0);
243

    
244
    /**
245
     * Returns a random variate from the lognormal distribution with "scale"
246
     * parameter m and "shape" parameter w. m and w correspond to the parameters
247
     * of the underlying normal distribution (m: mean, w: standard deviation.)
248
     *
249
     * Generation is using relationship to normal distribution.
250
     *
251
     * @param m  "scale" parameter, m>0
252
     * @param w  "shape" parameter, w>0
253
     * @param rng the underlying random number generator
254
     */
255
    inline double lognormal(double m, double w, int rng=0)
256
    {
257
        return exp(normal(m, w, rng));
258
    }
259

    
260
    /**
261
     * Returns a random variate from the Weibull distribution with parameters
262
     * a, b > 0, where a is the "scale" parameter and b is the "shape" parameter.
263
     * Sometimes Weibull is given with alpha and beta parameters, then alpha=b
264
     * and beta=a.
265
     *
266
     * The Weibull distribution gives the distribution of lifetimes of objects.
267
     * It was originally proposed to quantify fatigue data, but it is also used
268
     * in reliability analysis of systems involving a "weakest link," e.g.
269
     * in calculating a device's mean time to failure.
270
     *
271
     * When b=1, Weibull(a,b) is exponential with mean a.
272
     *
273
     * Generation uses inverse transform.
274
     *
275
     * @param a  the "scale" parameter, a>0
276
     * @param b  the "shape" parameter, b>0
277
     * @param rng the underlying random number generator
278
     */
279
    double weibull(double a, double b, int rng=0);
280

    
281
    /**
282
     * Returns a random variate from the shifted generalized Pareto distribution.
283
     *
284
     * Generation uses inverse transform.
285
     *
286
     * @param a,b  the usual parameters for generalized Pareto
287
     * @param c    shift parameter for left-shift
288
     * @param rng the underlying random number generator
289
     */
290
    double pareto_shifted(double a, double b, double c, int rng=0);
291

    
292
    //@}
293

    
294
    /**
295
     * @name Discrete distributions
296
     *
297
     * @ingroup RandomNumbers
298
     */
299
    //@{
300

    
301
    /**
302
     * Returns a random integer with uniform distribution in the range [a,b],
303
     * inclusive. (Note that the function can also return b.)
304
     *
305
     * @param a, b  the interval, a<b
306
     * @param rng the underlying random number generator
307
     */
308
    int intuniform(int a, int b, int rng=0);
309

    
310
    /**
311
     * Returns the result of a Bernoulli trial with probability p,
312
     * that is, 1 with probability p and 0 with probability (1-p).
313
     *
314
     * Generation is using elementary look-up.
315
     *
316
     * @param p  0=<p<=1
317
     * @param rng the underlying random number generator
318
     */
319
    inline int bernoulli(double p, int rng=0)
320
    {
321
        double U = genk_dblrand(rng);
322
        return (p > U) ? 1 : 0;
323
    }
324

    
325
    /**
326
     * Returns a random integer from the binomial distribution with
327
     * parameters n and p, that is, the number of successes in n independent
328
     * trials with probability p.
329
     *
330
     * Generation is using the relationship to Bernoulli distribution (runtime
331
     * is proportional to n).
332
     *
333
     * @param n n>=0
334
     * @param p 0<=p<=1
335
     * @param rng the underlying random number generator
336
     */
337
    int binomial(int n, double p, int rng=0);
338

    
339
    /**
340
     * Returns a random integer from the geometric distribution with parameter p,
341
     * that is, the number of independent trials with probability p until the
342
     * first success.
343
     *
344
     * This is the n=1 special case of the negative binomial distribution.
345
     *
346
     * Generation uses inverse transform.
347
     *
348
     * @param p  0<p<=1
349
     * @param rng the underlying random number generator
350
     */
351
    int geometric(double p, int rng=0);
352

    
353
    /**
354
     * Returns a random integer from the negative binomial distribution with
355
     * parameters n and p, that is, the number of failures occurring before
356
     * n successes in independent trials with probability p of success.
357
     *
358
     * Generation is using the relationship to geometric distribution (runtime is
359
     * proportional to n).
360
     *
361
     * @param n  n>=0
362
     * @param p  0<p<1
363
     * @param rng the underlying random number generator
364
     */
365
    int negbinomial(int n, double p, int rng=0);
366

    
367

    
368
    /**
369
     * Returns a random integer from the Poisson distribution with parameter lambda,
370
     * that is, the number of arrivals over unit time where the time between
371
     * successive arrivals follow exponential distribution with parameter lambda.
372
     *
373
     * Lambda is also the mean (and variance) of the distribution.
374
     *
375
     * Generation method depends on value of lambda:
376
     *
377
     *   - 0<lambda<=30: count number of events
378
     *   - lambda>30: Acceptance-Rejection due to Atkinson (see Banks, page 166)
379
     *
380
     * @param lambda  > 0
381
     * @param rng the underlying random number generator
382
     */
383
    int poisson(double lambda, int rng=0);
384

    
385
    /**
386
     * Produces random integer in range [0,r) using generator 0.
387
     */
388
    long intrand(long r);
389

    
390
    /**
391
     * Produces random double in range [0,1) using generator 0.
392
     */
393
    double dblrand();
394
    //@}
395

    
396
    /**
397
     * @name Compatibility
398
     *
399
     * @ingroup RandomNumbers
400
     */
401
    //@{
402
    /**
403
     * DEPRECATED: use uniform() instead.
404
     */
405
    double genk_uniform(double gen_nr, double a, double b);
406

    
407
    /**
408
     * DEPRECATED: use intuniform() instead.
409
     */
410
    double genk_intuniform(double gen_nr, double a, double b);
411

    
412
    /**
413
     * DEPRECATED: use exponential() instead.
414
     */
415
    double genk_exponential(double gen_nr, double p);
416

    
417
    /**
418
     * DEPRECATED: use normal() instead.
419
     */
420
    double genk_normal(double gen_nr, double mean, double variance);
421

    
422
    /**
423
     * DEPRECATED: use truncnormal() instead.
424
     */
425
    double genk_truncnormal(double gen_nr, double mean, double variance);
426
    //@}
427
};
428

    
429
#endif  // __CNUMGEN_H_