Branch data Line data Source code
1 : : // Copyright (C) 2010 Internet Systems Consortium, Inc. ("ISC")
2 : : //
3 : : // Permission to use, copy, modify, and/or distribute this software for any
4 : : // purpose with or without fee is hereby granted, provided that the above
5 : : // copyright notice and this permission notice appear in all copies.
6 : : //
7 : : // THE SOFTWARE IS PROVIDED "AS IS" AND ISC DISCLAIMS ALL WARRANTIES WITH
8 : : // REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED WARRANTIES OF MERCHANTABILITY
9 : : // AND FITNESS. IN NO EVENT SHALL ISC BE LIABLE FOR ANY SPECIAL, DIRECT,
10 : : // INDIRECT, OR CONSEQUENTIAL DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM
11 : : // LOSS OF USE, DATA OR PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE
12 : : // OR OTHER TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR
13 : : // PERFORMANCE OF THIS SOFTWARE.
14 : :
15 : : #ifndef __NSAS_RANDOM_NUMBER_GENERATOR_H
16 : : #define __NSAS_RANDOM_NUMBER_GENERATOR_H
17 : :
18 : : #include <algorithm>
19 : : #include <cmath>
20 : : #include <numeric>
21 : :
22 : : #include <exceptions/exceptions.h>
23 : :
24 : : #include <boost/random/mersenne_twister.hpp>
25 : : #include <boost/random/uniform_int.hpp>
26 : : #include <boost/random/uniform_real.hpp>
27 : : #include <boost/random/variate_generator.hpp>
28 : :
29 : : namespace isc {
30 : : namespace util {
31 : : namespace random {
32 : :
33 : 1 : class InvalidLimits : public isc::BadValue {
34 : : public:
35 : : InvalidLimits(const char* file, size_t line, const char* what) :
36 [ + - ]: 1 : isc::BadValue(file, line, what) {}
37 : : };
38 : :
39 : 2 : class SumNotOne : public isc::BadValue {
40 : : public:
41 : : SumNotOne(const char* file, size_t line, const char* what) :
42 [ + - ]: 2 : isc::BadValue(file, line, what) {}
43 : : };
44 : :
45 : 2 : class InvalidProbValue : public isc::BadValue {
46 : : public:
47 : : InvalidProbValue(const char* file, size_t line, const char* what) :
48 [ + - ]: 2 : isc::BadValue(file, line, what) {}
49 : : };
50 : :
51 : :
52 : :
53 : : /// \brief Uniform random integer generator
54 : : ///
55 : : /// Generate uniformly distributed integers in range of [min, max]
56 : : class UniformRandomIntegerGenerator{
57 : : public:
58 : : /// \brief Constructor
59 : : ///
60 : : /// \param min The minimum number in the range
61 : : /// \param max The maximum number in the range
62 : 3 : UniformRandomIntegerGenerator(int min, int max):
63 : : min_(std::min(min, max)), max_(std::max(min, max)),
64 : 12 : dist_(min_, max_), generator_(rng_, dist_)
65 : : {
66 : : // To preserve the restriction of the underlying uniform_int class (and
67 : : // to retain compatibility with earlier versions of the class), we will
68 : : // abort if the minimum and maximum given are the wrong way round.
69 [ + + ]: 3 : if (min > max) {
70 [ + - ]: 2 : isc_throw(InvalidLimits, "minimum limit is greater than maximum "
71 : : "when initializing UniformRandomIntegerGenerator");
72 : : }
73 : :
74 : : // Init with the current time
75 : 2 : rng_.seed(time(NULL));
76 : 2 : }
77 : :
78 : : /// \brief Generate uniformly distributed integer
79 : 100 : int operator()() { return generator_(); }
80 : : private:
81 : : /// Hide default and copy constructor
82 : : UniformRandomIntegerGenerator();///< Default constructor
83 : : UniformRandomIntegerGenerator(const UniformRandomIntegerGenerator&); ///< Copy constructor
84 : :
85 : : int min_; ///< The minimum integer that can generate
86 : : int max_; ///< The maximum integer that can generate
87 : : boost::uniform_int<> dist_; ///< Distribute uniformly.
88 : : boost::mt19937 rng_; ///< Mersenne Twister: A 623-dimensionally equidistributed uniform pseudo-random number generator
89 : : boost::variate_generator<boost::mt19937&, boost::uniform_int<> > generator_; ///< Uniform generator
90 : : };
91 : :
92 : : /// \brief Weighted random integer generator
93 : : ///
94 : : /// Generate random integers according different probabilities
95 : 96 : class WeightedRandomIntegerGenerator {
96 : : public:
97 : : /// \brief Constructor
98 : : ///
99 : : /// \param probabilities The probabies for all the integers, the probability must be
100 : : /// between 0 and 1.0, the sum of probabilities must be equal to 1.
101 : : /// For example, if the probabilities contains the following values:
102 : : /// 0.5 0.3 0.2, the 1st integer will be generated more frequently than the
103 : : /// other integers and the probability is proportional to its value.
104 : : /// \param min The minimum integer that generated, other integers will be
105 : : /// min, min + 1, ..., min + probabilities.size() - 1
106 : 10 : WeightedRandomIntegerGenerator(const std::vector<double>& probabilities,
107 : : size_t min = 0):
108 : 4 : dist_(0, 1.0), uniform_real_gen_(rng_, dist_), min_(min)
109 : : {
110 : : // The probabilities must be valid. Checking is quite an expensive
111 : : // operation, so is only done in a debug build.
112 [ + + ][ - + ]: 10 : assert(areProbabilitiesValid(probabilities));
113 : :
114 : : // Calculate the partial sum of probabilities
115 : : std::partial_sum(probabilities.begin(), probabilities.end(),
116 [ + - ]: 6 : std::back_inserter(cumulative_));
117 : : // Init with the current time
118 : 6 : rng_.seed(time(NULL));
119 : 6 : }
120 : :
121 : : /// \brief Default constructor
122 : : ///
123 : : WeightedRandomIntegerGenerator():
124 : 460 : dist_(0, 1.0), uniform_real_gen_(rng_, dist_), min_(0)
125 : : {
126 : : }
127 : :
128 : : /// \brief Reset the probabilities
129 : : ///
130 : : /// Change the weights of each integers
131 : : /// \param probabilities The probabies for all the integers
132 : : /// \param min The minimum integer that generated
133 : : void reset(const std::vector<double>& probabilities, size_t min = 0)
134 : : {
135 : : // The probabilities must be valid.
136 [ + - ][ - + ]: 402028 : assert(areProbabilitiesValid(probabilities));
137 : :
138 : : // Reset the cumulative sum
139 : 402028 : cumulative_.clear();
140 : :
141 : : // Calculate the partial sum of probabilities
142 : : std::partial_sum(probabilities.begin(), probabilities.end(),
143 [ + - ]: 402028 : std::back_inserter(cumulative_));
144 : :
145 : : // Reset the minimum integer
146 : 402028 : min_ = min;
147 : : }
148 : :
149 : : /// \brief Generate weighted random integer
150 : 600100 : size_t operator()()
151 : : {
152 : 3006387 : return std::lower_bound(cumulative_.begin(), cumulative_.end(), uniform_real_gen_())
153 : 1806187 : - cumulative_.begin() + min_;
154 : : }
155 : :
156 : : private:
157 : : /// \brief Check the validation of probabilities vector
158 : : ///
159 : : /// The probability must be in range of [0, 1.0] and the sum must be equal
160 : : /// to 1.0. Empty probabilities are also valid.
161 : : ///
162 : : /// Checking the probabilities is quite an expensive operation, so it is
163 : : /// only done during a debug build (via a call through assert()). However,
164 : : /// instead of letting assert() call abort(), if this method encounters an
165 : : /// error, an exception is thrown. This makes unit testing somewhat easier.
166 : : ///
167 : : /// \param probabilities Vector of probabilities.
168 : 402038 : bool areProbabilitiesValid(const std::vector<double>& probabilities) const
169 : : {
170 : : typedef std::vector<double>::const_iterator Iterator;
171 [ + + ]: 402038 : double sum = probabilities.empty() ? 1.0 : 0.0;
172 [ + + ]: 1606092 : for(Iterator it = probabilities.begin(); it != probabilities.end(); ++it){
173 : : //The probability must be in [0, 1.0]
174 [ + + ][ + + ]: 1204056 : if(*it < 0.0 || *it > 1.0) {
[ + + ]
175 [ + - ]: 4 : isc_throw(InvalidProbValue,
176 : : "probability must be in the range 0..1");
177 : : }
178 : :
179 : 1204054 : sum += *it;
180 : : }
181 : :
182 : 402036 : double epsilon = 0.0001;
183 : : // The sum must be equal to 1
184 [ + + ]: 402036 : if (std::fabs(sum - 1.0) >= epsilon) {
185 [ + - ]: 4 : isc_throw(SumNotOne, "Sum of probabilities is not equal to 1");
186 : : }
187 : :
188 : 402034 : return true;
189 : : }
190 : :
191 : : std::vector<double> cumulative_; ///< Partial sum of the probabilities
192 : : boost::mt19937 rng_; ///< Mersenne Twister: A 623-dimensionally equidistributed uniform pseudo-random number generator
193 : : boost::uniform_real<> dist_; ///< Uniformly distributed real numbers
194 : :
195 : : // Shortcut typedef
196 : : // This typedef is placed directly before its use, as the sunstudio
197 : : // compiler could not handle it being anywhere else (don't know why)
198 : : typedef boost::variate_generator<boost::mt19937&, boost::uniform_real<> > UniformRealGenerator;
199 : : UniformRealGenerator uniform_real_gen_; ///< Uniformly distributed random real numbers generator
200 : :
201 : : size_t min_; ///< The minimum integer that will be generated
202 : : };
203 : :
204 : : } // namespace random
205 : : } // namespace util
206 : : } // namespace isc
207 : :
208 : : #endif//__NSAS_RANDOM_NUMBER_GENERATOR_H
|
