1 //===----------------------------------------------------------------------===//
2 //
3 // Part of the LLVM Project, under the Apache License v2.0 with LLVM Exceptions.
4 // See https://llvm.org/LICENSE.txt for license information.
5 // SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception
6 //
7 //===----------------------------------------------------------------------===//
8 //
9 // REQUIRES: long_tests
10 
11 // <random>
12 
13 // template<class IntType = int>
14 // class geometric_distribution
15 
16 // template<class _URNG> result_type operator()(_URNG& g);
17 
18 #include <random>
19 #include <numeric>
20 #include <vector>
21 #include <cassert>
22 
23 #include "test_macros.h"
24 
25 template <class T>
26 inline
27 T
sqr(T x)28 sqr(T x)
29 {
30     return x * x;
31 }
32 
test_small_inputs()33 void test_small_inputs() {
34   std::mt19937 engine;
35   std::geometric_distribution<std::int16_t> distribution(5.45361e-311);
36   typedef std::geometric_distribution<std::int16_t>::result_type result_type;
37   for (int i = 0; i < 1000; ++i) {
38     volatile result_type res = distribution(engine);
39     ((void)res);
40   }
41 }
42 
43 void
test1()44 test1()
45 {
46     typedef std::geometric_distribution<> D;
47     typedef std::mt19937 G;
48     G g;
49     D d(.03125);
50     const int N = 1000000;
51     std::vector<D::result_type> u;
52     for (int i = 0; i < N; ++i)
53     {
54         D::result_type v = d(g);
55         assert(d.min() <= v && v <= d.max());
56         u.push_back(v);
57     }
58     double mean = std::accumulate(u.begin(), u.end(),
59                                           double(0)) / u.size();
60     double var = 0;
61     double skew = 0;
62     double kurtosis = 0;
63     for (unsigned i = 0; i < u.size(); ++i)
64     {
65         double dbl = (u[i] - mean);
66         double d2 = sqr(dbl);
67         var += d2;
68         skew += dbl * d2;
69         kurtosis += d2 * d2;
70     }
71     var /= u.size();
72     double dev = std::sqrt(var);
73     skew /= u.size() * dev * var;
74     kurtosis /= u.size() * var * var;
75     kurtosis -= 3;
76     double x_mean = (1 - d.p()) / d.p();
77     double x_var = x_mean / d.p();
78     double x_skew = (2 - d.p()) / std::sqrt((1 - d.p()));
79     double x_kurtosis = 6 + sqr(d.p()) / (1 - d.p());
80     assert(std::abs((mean - x_mean) / x_mean) < 0.01);
81     assert(std::abs((var - x_var) / x_var) < 0.01);
82     assert(std::abs((skew - x_skew) / x_skew) < 0.01);
83     assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.01);
84 }
85 
86 void
test2()87 test2()
88 {
89     typedef std::geometric_distribution<> D;
90     typedef std::mt19937 G;
91     G g;
92     D d(0.05);
93     const int N = 1000000;
94     std::vector<D::result_type> u;
95     for (int i = 0; i < N; ++i)
96     {
97         D::result_type v = d(g);
98         assert(d.min() <= v && v <= d.max());
99         u.push_back(v);
100     }
101     double mean = std::accumulate(u.begin(), u.end(),
102                                           double(0)) / u.size();
103     double var = 0;
104     double skew = 0;
105     double kurtosis = 0;
106     for (unsigned i = 0; i < u.size(); ++i)
107     {
108         double dbl = (u[i] - mean);
109         double d2 = sqr(dbl);
110         var += d2;
111         skew += dbl * d2;
112         kurtosis += d2 * d2;
113     }
114     var /= u.size();
115     double dev = std::sqrt(var);
116     skew /= u.size() * dev * var;
117     kurtosis /= u.size() * var * var;
118     kurtosis -= 3;
119     double x_mean = (1 - d.p()) / d.p();
120     double x_var = x_mean / d.p();
121     double x_skew = (2 - d.p()) / std::sqrt((1 - d.p()));
122     double x_kurtosis = 6 + sqr(d.p()) / (1 - d.p());
123     assert(std::abs((mean - x_mean) / x_mean) < 0.01);
124     assert(std::abs((var - x_var) / x_var) < 0.01);
125     assert(std::abs((skew - x_skew) / x_skew) < 0.01);
126     assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.03);
127 }
128 
129 void
test3()130 test3()
131 {
132     typedef std::geometric_distribution<> D;
133     typedef std::minstd_rand G;
134     G g;
135     D d(.25);
136     const int N = 1000000;
137     std::vector<D::result_type> u;
138     for (int i = 0; i < N; ++i)
139     {
140         D::result_type v = d(g);
141         assert(d.min() <= v && v <= d.max());
142         u.push_back(v);
143     }
144     double mean = std::accumulate(u.begin(), u.end(),
145                                           double(0)) / u.size();
146     double var = 0;
147     double skew = 0;
148     double kurtosis = 0;
149     for (unsigned i = 0; i < u.size(); ++i)
150     {
151         double dbl = (u[i] - mean);
152         double d2 = sqr(dbl);
153         var += d2;
154         skew += dbl * d2;
155         kurtosis += d2 * d2;
156     }
157     var /= u.size();
158     double dev = std::sqrt(var);
159     skew /= u.size() * dev * var;
160     kurtosis /= u.size() * var * var;
161     kurtosis -= 3;
162     double x_mean = (1 - d.p()) / d.p();
163     double x_var = x_mean / d.p();
164     double x_skew = (2 - d.p()) / std::sqrt((1 - d.p()));
165     double x_kurtosis = 6 + sqr(d.p()) / (1 - d.p());
166     assert(std::abs((mean - x_mean) / x_mean) < 0.01);
167     assert(std::abs((var - x_var) / x_var) < 0.01);
168     assert(std::abs((skew - x_skew) / x_skew) < 0.01);
169     assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.02);
170 }
171 
172 void
test4()173 test4()
174 {
175     typedef std::geometric_distribution<> D;
176     typedef std::mt19937 G;
177     G g;
178     D d(0.5);
179     const int N = 1000000;
180     std::vector<D::result_type> u;
181     for (int i = 0; i < N; ++i)
182     {
183         D::result_type v = d(g);
184         assert(d.min() <= v && v <= d.max());
185         u.push_back(v);
186     }
187     double mean = std::accumulate(u.begin(), u.end(),
188                                           double(0)) / u.size();
189     double var = 0;
190     double skew = 0;
191     double kurtosis = 0;
192     for (unsigned i = 0; i < u.size(); ++i)
193     {
194         double dbl = (u[i] - mean);
195         double d2 = sqr(dbl);
196         var += d2;
197         skew += dbl * d2;
198         kurtosis += d2 * d2;
199     }
200     var /= u.size();
201     double dev = std::sqrt(var);
202     skew /= u.size() * dev * var;
203     kurtosis /= u.size() * var * var;
204     kurtosis -= 3;
205     double x_mean = (1 - d.p()) / d.p();
206     double x_var = x_mean / d.p();
207     double x_skew = (2 - d.p()) / std::sqrt((1 - d.p()));
208     double x_kurtosis = 6 + sqr(d.p()) / (1 - d.p());
209     assert(std::abs((mean - x_mean) / x_mean) < 0.01);
210     assert(std::abs((var - x_var) / x_var) < 0.01);
211     assert(std::abs((skew - x_skew) / x_skew) < 0.01);
212     assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.02);
213 }
214 
215 void
test5()216 test5()
217 {
218     typedef std::geometric_distribution<> D;
219     typedef std::mt19937 G;
220     G g;
221     D d(0.75);
222     const int N = 1000000;
223     std::vector<D::result_type> u;
224     for (int i = 0; i < N; ++i)
225     {
226         D::result_type v = d(g);
227         assert(d.min() <= v && v <= d.max());
228         u.push_back(v);
229     }
230     double mean = std::accumulate(u.begin(), u.end(),
231                                           double(0)) / u.size();
232     double var = 0;
233     double skew = 0;
234     double kurtosis = 0;
235     for (unsigned i = 0; i < u.size(); ++i)
236     {
237         double dbl = (u[i] - mean);
238         double d2 = sqr(dbl);
239         var += d2;
240         skew += dbl * d2;
241         kurtosis += d2 * d2;
242     }
243     var /= u.size();
244     double dev = std::sqrt(var);
245     skew /= u.size() * dev * var;
246     kurtosis /= u.size() * var * var;
247     kurtosis -= 3;
248     double x_mean = (1 - d.p()) / d.p();
249     double x_var = x_mean / d.p();
250     double x_skew = (2 - d.p()) / std::sqrt((1 - d.p()));
251     double x_kurtosis = 6 + sqr(d.p()) / (1 - d.p());
252     assert(std::abs((mean - x_mean) / x_mean) < 0.01);
253     assert(std::abs((var - x_var) / x_var) < 0.01);
254     assert(std::abs((skew - x_skew) / x_skew) < 0.01);
255     assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.02);
256 }
257 
258 void
test6()259 test6()
260 {
261     typedef std::geometric_distribution<> D;
262     typedef std::mt19937 G;
263     G g;
264     D d(0.96875);
265     const int N = 1000000;
266     std::vector<D::result_type> u;
267     for (int i = 0; i < N; ++i)
268     {
269         D::result_type v = d(g);
270         assert(d.min() <= v && v <= d.max());
271         u.push_back(v);
272     }
273     double mean = std::accumulate(u.begin(), u.end(),
274                                           double(0)) / u.size();
275     double var = 0;
276     double skew = 0;
277     double kurtosis = 0;
278     for (unsigned i = 0; i < u.size(); ++i)
279     {
280         double dbl = (u[i] - mean);
281         double d2 = sqr(dbl);
282         var += d2;
283         skew += dbl * d2;
284         kurtosis += d2 * d2;
285     }
286     var /= u.size();
287     double dev = std::sqrt(var);
288     skew /= u.size() * dev * var;
289     kurtosis /= u.size() * var * var;
290     kurtosis -= 3;
291     double x_mean = (1 - d.p()) / d.p();
292     double x_var = x_mean / d.p();
293     double x_skew = (2 - d.p()) / std::sqrt((1 - d.p()));
294     double x_kurtosis = 6 + sqr(d.p()) / (1 - d.p());
295     assert(std::abs((mean - x_mean) / x_mean) < 0.01);
296     assert(std::abs((var - x_var) / x_var) < 0.01);
297     assert(std::abs((skew - x_skew) / x_skew) < 0.01);
298     assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.02);
299 }
300 
main(int,char **)301 int main(int, char**)
302 {
303     test1();
304     test2();
305     test3();
306     test4();
307     test5();
308     test6();
309     test_small_inputs();
310 
311   return 0;
312 }
313