1 /*
2  *  Copyright (c) 2015 The WebRTC project authors. All Rights Reserved.
3  *
4  *  Use of this source code is governed by a BSD-style license
5  *  that can be found in the LICENSE file in the root of the source
6  *  tree. An additional intellectual property rights grant can be found
7  *  in the file PATENTS.  All contributing project authors may
8  *  be found in the AUTHORS file in the root of the source tree.
9  */
10 
11 #include "rtc_base/random.h"
12 
13 #include <math.h>
14 
15 #include <limits>
16 #include <vector>
17 
18 #include "rtc_base/numerics/math_utils.h"  // unsigned difference
19 #include "test/gtest.h"
20 
21 namespace webrtc {
22 
23 namespace {
24 // Computes the positive remainder of x/n.
25 template <typename T>
fdiv_remainder(T x,T n)26 T fdiv_remainder(T x, T n) {
27   RTC_CHECK_GE(n, 0);
28   T remainder = x % n;
29   if (remainder < 0)
30     remainder += n;
31   return remainder;
32 }
33 }  // namespace
34 
35 // Sample a number of random integers of type T. Divide them into buckets
36 // based on the remainder when dividing by bucket_count and check that each
37 // bucket gets roughly the expected number of elements.
38 template <typename T>
UniformBucketTest(T bucket_count,int samples,Random * prng)39 void UniformBucketTest(T bucket_count, int samples, Random* prng) {
40   std::vector<int> buckets(bucket_count, 0);
41 
42   uint64_t total_values = 1ull << (std::numeric_limits<T>::digits +
43                                    std::numeric_limits<T>::is_signed);
44   T upper_limit =
45       std::numeric_limits<T>::max() -
46       static_cast<T>(total_values % static_cast<uint64_t>(bucket_count));
47   ASSERT_GT(upper_limit, std::numeric_limits<T>::max() / 2);
48 
49   for (int i = 0; i < samples; i++) {
50     T sample;
51     do {
52       // We exclude a few numbers from the range so that it is divisible by
53       // the number of buckets. If we are unlucky and hit one of the excluded
54       // numbers we just resample. Note that if the number of buckets is a
55       // power of 2, then we don't have to exclude anything.
56       sample = prng->Rand<T>();
57     } while (sample > upper_limit);
58     buckets[fdiv_remainder(sample, bucket_count)]++;
59   }
60 
61   for (T i = 0; i < bucket_count; i++) {
62     // Expect the result to be within 3 standard deviations of the mean.
63     EXPECT_NEAR(buckets[i], samples / bucket_count,
64                 3 * sqrt(samples / bucket_count));
65   }
66 }
67 
TEST(RandomNumberGeneratorTest,BucketTestSignedChar)68 TEST(RandomNumberGeneratorTest, BucketTestSignedChar) {
69   Random prng(7297352569824ull);
70   UniformBucketTest<signed char>(64, 640000, &prng);
71   UniformBucketTest<signed char>(11, 440000, &prng);
72   UniformBucketTest<signed char>(3, 270000, &prng);
73 }
74 
TEST(RandomNumberGeneratorTest,BucketTestUnsignedChar)75 TEST(RandomNumberGeneratorTest, BucketTestUnsignedChar) {
76   Random prng(7297352569824ull);
77   UniformBucketTest<unsigned char>(64, 640000, &prng);
78   UniformBucketTest<unsigned char>(11, 440000, &prng);
79   UniformBucketTest<unsigned char>(3, 270000, &prng);
80 }
81 
TEST(RandomNumberGeneratorTest,BucketTestSignedShort)82 TEST(RandomNumberGeneratorTest, BucketTestSignedShort) {
83   Random prng(7297352569824ull);
84   UniformBucketTest<int16_t>(64, 640000, &prng);
85   UniformBucketTest<int16_t>(11, 440000, &prng);
86   UniformBucketTest<int16_t>(3, 270000, &prng);
87 }
88 
TEST(RandomNumberGeneratorTest,BucketTestUnsignedShort)89 TEST(RandomNumberGeneratorTest, BucketTestUnsignedShort) {
90   Random prng(7297352569824ull);
91   UniformBucketTest<uint16_t>(64, 640000, &prng);
92   UniformBucketTest<uint16_t>(11, 440000, &prng);
93   UniformBucketTest<uint16_t>(3, 270000, &prng);
94 }
95 
TEST(RandomNumberGeneratorTest,BucketTestSignedInt)96 TEST(RandomNumberGeneratorTest, BucketTestSignedInt) {
97   Random prng(7297352569824ull);
98   UniformBucketTest<signed int>(64, 640000, &prng);
99   UniformBucketTest<signed int>(11, 440000, &prng);
100   UniformBucketTest<signed int>(3, 270000, &prng);
101 }
102 
TEST(RandomNumberGeneratorTest,BucketTestUnsignedInt)103 TEST(RandomNumberGeneratorTest, BucketTestUnsignedInt) {
104   Random prng(7297352569824ull);
105   UniformBucketTest<unsigned int>(64, 640000, &prng);
106   UniformBucketTest<unsigned int>(11, 440000, &prng);
107   UniformBucketTest<unsigned int>(3, 270000, &prng);
108 }
109 
110 // The range of the random numbers is divided into bucket_count intervals
111 // of consecutive numbers. Check that approximately equally many numbers
112 // from each inteval are generated.
BucketTestSignedInterval(unsigned int bucket_count,unsigned int samples,int32_t low,int32_t high,int sigma_level,Random * prng)113 void BucketTestSignedInterval(unsigned int bucket_count,
114                               unsigned int samples,
115                               int32_t low,
116                               int32_t high,
117                               int sigma_level,
118                               Random* prng) {
119   std::vector<unsigned int> buckets(bucket_count, 0);
120 
121   ASSERT_GE(high, low);
122   ASSERT_GE(bucket_count, 2u);
123   uint32_t interval = unsigned_difference<int32_t>(high, low) + 1;
124   uint32_t numbers_per_bucket;
125   if (interval == 0) {
126     // The computation high - low + 1 should be 2^32 but overflowed
127     // Hence, bucket_count must be a power of 2
128     ASSERT_EQ(bucket_count & (bucket_count - 1), 0u);
129     numbers_per_bucket = (0x80000000u / bucket_count) * 2;
130   } else {
131     ASSERT_EQ(interval % bucket_count, 0u);
132     numbers_per_bucket = interval / bucket_count;
133   }
134 
135   for (unsigned int i = 0; i < samples; i++) {
136     int32_t sample = prng->Rand(low, high);
137     EXPECT_LE(low, sample);
138     EXPECT_GE(high, sample);
139     buckets[unsigned_difference<int32_t>(sample, low) / numbers_per_bucket]++;
140   }
141 
142   for (unsigned int i = 0; i < bucket_count; i++) {
143     // Expect the result to be within 3 standard deviations of the mean,
144     // or more generally, within sigma_level standard deviations of the mean.
145     double mean = static_cast<double>(samples) / bucket_count;
146     EXPECT_NEAR(buckets[i], mean, sigma_level * sqrt(mean));
147   }
148 }
149 
150 // The range of the random numbers is divided into bucket_count intervals
151 // of consecutive numbers. Check that approximately equally many numbers
152 // from each inteval are generated.
BucketTestUnsignedInterval(unsigned int bucket_count,unsigned int samples,uint32_t low,uint32_t high,int sigma_level,Random * prng)153 void BucketTestUnsignedInterval(unsigned int bucket_count,
154                                 unsigned int samples,
155                                 uint32_t low,
156                                 uint32_t high,
157                                 int sigma_level,
158                                 Random* prng) {
159   std::vector<unsigned int> buckets(bucket_count, 0);
160 
161   ASSERT_GE(high, low);
162   ASSERT_GE(bucket_count, 2u);
163   uint32_t interval = high - low + 1;
164   uint32_t numbers_per_bucket;
165   if (interval == 0) {
166     // The computation high - low + 1 should be 2^32 but overflowed
167     // Hence, bucket_count must be a power of 2
168     ASSERT_EQ(bucket_count & (bucket_count - 1), 0u);
169     numbers_per_bucket = (0x80000000u / bucket_count) * 2;
170   } else {
171     ASSERT_EQ(interval % bucket_count, 0u);
172     numbers_per_bucket = interval / bucket_count;
173   }
174 
175   for (unsigned int i = 0; i < samples; i++) {
176     uint32_t sample = prng->Rand(low, high);
177     EXPECT_LE(low, sample);
178     EXPECT_GE(high, sample);
179     buckets[(sample - low) / numbers_per_bucket]++;
180   }
181 
182   for (unsigned int i = 0; i < bucket_count; i++) {
183     // Expect the result to be within 3 standard deviations of the mean,
184     // or more generally, within sigma_level standard deviations of the mean.
185     double mean = static_cast<double>(samples) / bucket_count;
186     EXPECT_NEAR(buckets[i], mean, sigma_level * sqrt(mean));
187   }
188 }
189 
TEST(RandomNumberGeneratorTest,UniformUnsignedInterval)190 TEST(RandomNumberGeneratorTest, UniformUnsignedInterval) {
191   Random prng(299792458ull);
192   BucketTestUnsignedInterval(2, 100000, 0, 1, 3, &prng);
193   BucketTestUnsignedInterval(7, 100000, 1, 14, 3, &prng);
194   BucketTestUnsignedInterval(11, 100000, 1000, 1010, 3, &prng);
195   BucketTestUnsignedInterval(100, 100000, 0, 99, 3, &prng);
196   BucketTestUnsignedInterval(2, 100000, 0, 4294967295, 3, &prng);
197   BucketTestUnsignedInterval(17, 100000, 455, 2147484110, 3, &prng);
198   // 99.7% of all samples will be within 3 standard deviations of the mean,
199   // but since we test 1000 buckets we allow an interval of 4 sigma.
200   BucketTestUnsignedInterval(1000, 1000000, 0, 2147483999, 4, &prng);
201 }
202 
TEST(RandomNumberGeneratorTest,UniformSignedInterval)203 TEST(RandomNumberGeneratorTest, UniformSignedInterval) {
204   Random prng(66260695729ull);
205   BucketTestSignedInterval(2, 100000, 0, 1, 3, &prng);
206   BucketTestSignedInterval(7, 100000, -2, 4, 3, &prng);
207   BucketTestSignedInterval(11, 100000, 1000, 1010, 3, &prng);
208   BucketTestSignedInterval(100, 100000, 0, 99, 3, &prng);
209   BucketTestSignedInterval(2, 100000, std::numeric_limits<int32_t>::min(),
210                            std::numeric_limits<int32_t>::max(), 3, &prng);
211   BucketTestSignedInterval(17, 100000, -1073741826, 1073741829, 3, &prng);
212   // 99.7% of all samples will be within 3 standard deviations of the mean,
213   // but since we test 1000 buckets we allow an interval of 4 sigma.
214   BucketTestSignedInterval(1000, 1000000, -352, 2147483647, 4, &prng);
215 }
216 
217 // The range of the random numbers is divided into bucket_count intervals
218 // of consecutive numbers. Check that approximately equally many numbers
219 // from each inteval are generated.
BucketTestFloat(unsigned int bucket_count,unsigned int samples,int sigma_level,Random * prng)220 void BucketTestFloat(unsigned int bucket_count,
221                      unsigned int samples,
222                      int sigma_level,
223                      Random* prng) {
224   ASSERT_GE(bucket_count, 2u);
225   std::vector<unsigned int> buckets(bucket_count, 0);
226 
227   for (unsigned int i = 0; i < samples; i++) {
228     uint32_t sample = bucket_count * prng->Rand<float>();
229     EXPECT_LE(0u, sample);
230     EXPECT_GE(bucket_count - 1, sample);
231     buckets[sample]++;
232   }
233 
234   for (unsigned int i = 0; i < bucket_count; i++) {
235     // Expect the result to be within 3 standard deviations of the mean,
236     // or more generally, within sigma_level standard deviations of the mean.
237     double mean = static_cast<double>(samples) / bucket_count;
238     EXPECT_NEAR(buckets[i], mean, sigma_level * sqrt(mean));
239   }
240 }
241 
TEST(RandomNumberGeneratorTest,UniformFloatInterval)242 TEST(RandomNumberGeneratorTest, UniformFloatInterval) {
243   Random prng(1380648813ull);
244   BucketTestFloat(100, 100000, 3, &prng);
245   // 99.7% of all samples will be within 3 standard deviations of the mean,
246   // but since we test 1000 buckets we allow an interval of 4 sigma.
247   // BucketTestSignedInterval(1000, 1000000, -352, 2147483647, 4, &prng);
248 }
249 
TEST(RandomNumberGeneratorTest,SignedHasSameBitPattern)250 TEST(RandomNumberGeneratorTest, SignedHasSameBitPattern) {
251   Random prng_signed(66738480ull), prng_unsigned(66738480ull);
252 
253   for (int i = 0; i < 1000; i++) {
254     signed int s = prng_signed.Rand<signed int>();
255     unsigned int u = prng_unsigned.Rand<unsigned int>();
256     EXPECT_EQ(u, static_cast<unsigned int>(s));
257   }
258 
259   for (int i = 0; i < 1000; i++) {
260     int16_t s = prng_signed.Rand<int16_t>();
261     uint16_t u = prng_unsigned.Rand<uint16_t>();
262     EXPECT_EQ(u, static_cast<uint16_t>(s));
263   }
264 
265   for (int i = 0; i < 1000; i++) {
266     signed char s = prng_signed.Rand<signed char>();
267     unsigned char u = prng_unsigned.Rand<unsigned char>();
268     EXPECT_EQ(u, static_cast<unsigned char>(s));
269   }
270 }
271 
TEST(RandomNumberGeneratorTest,Gaussian)272 TEST(RandomNumberGeneratorTest, Gaussian) {
273   const int kN = 100000;
274   const int kBuckets = 100;
275   const double kMean = 49;
276   const double kStddev = 10;
277 
278   Random prng(1256637061);
279 
280   std::vector<unsigned int> buckets(kBuckets, 0);
281   for (int i = 0; i < kN; i++) {
282     int index = prng.Gaussian(kMean, kStddev) + 0.5;
283     if (index >= 0 && index < kBuckets) {
284       buckets[index]++;
285     }
286   }
287 
288   const double kPi = 3.14159265358979323846;
289   const double kScale = 1 / (kStddev * sqrt(2.0 * kPi));
290   const double kDiv = -2.0 * kStddev * kStddev;
291   for (int n = 0; n < kBuckets; ++n) {
292     // Use Simpsons rule to estimate the probability that a random gaussian
293     // sample is in the interval [n-0.5, n+0.5].
294     double f_left = kScale * exp((n - kMean - 0.5) * (n - kMean - 0.5) / kDiv);
295     double f_mid = kScale * exp((n - kMean) * (n - kMean) / kDiv);
296     double f_right = kScale * exp((n - kMean + 0.5) * (n - kMean + 0.5) / kDiv);
297     double normal_dist = (f_left + 4 * f_mid + f_right) / 6;
298     // Expect the number of samples to be within 3 standard deviations
299     // (rounded up) of the expected number of samples in the bucket.
300     EXPECT_NEAR(buckets[n], kN * normal_dist, 3 * sqrt(kN * normal_dist) + 1);
301   }
302 }
303 
304 }  // namespace webrtc
305