xref: /llvm-project/libcxx/test/std/numerics/rand/rand.dist/rand.dist.uni/rand.dist.uni.real/eval.pass.cpp (revision 09e3a360581dc36d0820d3fb6da9bd7cfed87b5d)
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 RealType = double>
14 // class uniform_real_distribution
15 
16 // template<class _URNG> result_type operator()(_URNG& g);
17 
18 #include <random>
19 #include <cassert>
20 #include <cmath>
21 #include <cstddef>
22 #include <numeric>
23 #include <vector>
24 
25 #include "test_macros.h"
26 
27 template <class T>
28 inline
29 T
30 sqr(T x)
31 {
32     return x * x;
33 }
34 
35 int main(int, char**)
36 {
37     {
38         typedef std::uniform_real_distribution<> D;
39         typedef std::minstd_rand0 G;
40         G g;
41         D d;
42         const int N = 100000;
43         std::vector<D::result_type> u;
44         for (int i = 0; i < N; ++i)
45         {
46             D::result_type v = d(g);
47             assert(d.a() <= v && v < d.b());
48             u.push_back(v);
49         }
50         D::result_type mean = std::accumulate(u.begin(), u.end(),
51                                               D::result_type(0)) / u.size();
52         D::result_type var = 0;
53         D::result_type skew = 0;
54         D::result_type kurtosis = 0;
55         for (std::size_t i = 0; i < u.size(); ++i)
56         {
57             D::result_type dbl = (u[i] - mean);
58             D::result_type d2 = sqr(dbl);
59             var += d2;
60             skew += dbl * d2;
61             kurtosis += d2 * d2;
62         }
63         var /= u.size();
64         D::result_type dev = std::sqrt(var);
65         skew /= u.size() * dev * var;
66         kurtosis /= u.size() * var * var;
67         kurtosis -= 3;
68         D::result_type x_mean = (d.a() + d.b()) / 2;
69         D::result_type x_var = sqr(d.b() - d.a()) / 12;
70         D::result_type x_skew = 0;
71         D::result_type x_kurtosis = -6./5;
72         assert(std::abs((mean - x_mean) / x_mean) < 0.01);
73         assert(std::abs((var - x_var) / x_var) < 0.01);
74         assert(std::abs(skew - x_skew) < 0.01);
75         assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.01);
76     }
77     {
78         typedef std::uniform_real_distribution<> D;
79         typedef std::minstd_rand G;
80         G g;
81         D d;
82         const int N = 100000;
83         std::vector<D::result_type> u;
84         for (int i = 0; i < N; ++i)
85         {
86             D::result_type v = d(g);
87             assert(d.a() <= v && v < d.b());
88             u.push_back(v);
89         }
90         D::result_type mean = std::accumulate(u.begin(), u.end(),
91                                               D::result_type(0)) / u.size();
92         D::result_type var = 0;
93         D::result_type skew = 0;
94         D::result_type kurtosis = 0;
95         for (std::size_t i = 0; i < u.size(); ++i)
96         {
97             D::result_type dbl = (u[i] - mean);
98             D::result_type d2 = sqr(dbl);
99             var += d2;
100             skew += dbl * d2;
101             kurtosis += d2 * d2;
102         }
103         var /= u.size();
104         D::result_type dev = std::sqrt(var);
105         skew /= u.size() * dev * var;
106         kurtosis /= u.size() * var * var;
107         kurtosis -= 3;
108         D::result_type x_mean = (d.a() + d.b()) / 2;
109         D::result_type x_var = sqr(d.b() - d.a()) / 12;
110         D::result_type x_skew = 0;
111         D::result_type x_kurtosis = -6./5;
112         assert(std::abs((mean - x_mean) / x_mean) < 0.01);
113         assert(std::abs((var - x_var) / x_var) < 0.01);
114         assert(std::abs(skew - x_skew) < 0.01);
115         assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.01);
116     }
117     {
118         typedef std::uniform_real_distribution<> D;
119         typedef std::mt19937 G;
120         G g;
121         D d;
122         const int N = 100000;
123         std::vector<D::result_type> u;
124         for (int i = 0; i < N; ++i)
125         {
126             D::result_type v = d(g);
127             assert(d.a() <= v && v < d.b());
128             u.push_back(v);
129         }
130         D::result_type mean = std::accumulate(u.begin(), u.end(),
131                                               D::result_type(0)) / u.size();
132         D::result_type var = 0;
133         D::result_type skew = 0;
134         D::result_type kurtosis = 0;
135         for (std::size_t i = 0; i < u.size(); ++i)
136         {
137             D::result_type dbl = (u[i] - mean);
138             D::result_type d2 = sqr(dbl);
139             var += d2;
140             skew += dbl * d2;
141             kurtosis += d2 * d2;
142         }
143         var /= u.size();
144         D::result_type dev = std::sqrt(var);
145         skew /= u.size() * dev * var;
146         kurtosis /= u.size() * var * var;
147         kurtosis -= 3;
148         D::result_type x_mean = (d.a() + d.b()) / 2;
149         D::result_type x_var = sqr(d.b() - d.a()) / 12;
150         D::result_type x_skew = 0;
151         D::result_type x_kurtosis = -6./5;
152         assert(std::abs((mean - x_mean) / x_mean) < 0.01);
153         assert(std::abs((var - x_var) / x_var) < 0.01);
154         assert(std::abs(skew - x_skew) < 0.01);
155         assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.01);
156     }
157     {
158         typedef std::uniform_real_distribution<> D;
159         typedef std::mt19937_64 G;
160         G g;
161         D d;
162         const int N = 100000;
163         std::vector<D::result_type> u;
164         for (int i = 0; i < N; ++i)
165         {
166             D::result_type v = d(g);
167             assert(d.a() <= v && v < d.b());
168             u.push_back(v);
169         }
170         D::result_type mean = std::accumulate(u.begin(), u.end(),
171                                               D::result_type(0)) / u.size();
172         D::result_type var = 0;
173         D::result_type skew = 0;
174         D::result_type kurtosis = 0;
175         for (std::size_t i = 0; i < u.size(); ++i)
176         {
177             D::result_type dbl = (u[i] - mean);
178             D::result_type d2 = sqr(dbl);
179             var += d2;
180             skew += dbl * d2;
181             kurtosis += d2 * d2;
182         }
183         var /= u.size();
184         D::result_type dev = std::sqrt(var);
185         skew /= u.size() * dev * var;
186         kurtosis /= u.size() * var * var;
187         kurtosis -= 3;
188         D::result_type x_mean = (d.a() + d.b()) / 2;
189         D::result_type x_var = sqr(d.b() - d.a()) / 12;
190         D::result_type x_skew = 0;
191         D::result_type x_kurtosis = -6./5;
192         assert(std::abs((mean - x_mean) / x_mean) < 0.01);
193         assert(std::abs((var - x_var) / x_var) < 0.01);
194         assert(std::abs(skew - x_skew) < 0.01);
195         assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.01);
196     }
197     {
198         typedef std::uniform_real_distribution<> D;
199         typedef std::ranlux24_base G;
200         G g;
201         D d;
202         const int N = 100000;
203         std::vector<D::result_type> u;
204         for (int i = 0; i < N; ++i)
205         {
206             D::result_type v = d(g);
207             assert(d.a() <= v && v < d.b());
208             u.push_back(v);
209         }
210         D::result_type mean = std::accumulate(u.begin(), u.end(),
211                                               D::result_type(0)) / u.size();
212         D::result_type var = 0;
213         D::result_type skew = 0;
214         D::result_type kurtosis = 0;
215         for (std::size_t i = 0; i < u.size(); ++i)
216         {
217             D::result_type dbl = (u[i] - mean);
218             D::result_type d2 = sqr(dbl);
219             var += d2;
220             skew += dbl * d2;
221             kurtosis += d2 * d2;
222         }
223         var /= u.size();
224         D::result_type dev = std::sqrt(var);
225         skew /= u.size() * dev * var;
226         kurtosis /= u.size() * var * var;
227         kurtosis -= 3;
228         D::result_type x_mean = (d.a() + d.b()) / 2;
229         D::result_type x_var = sqr(d.b() - d.a()) / 12;
230         D::result_type x_skew = 0;
231         D::result_type x_kurtosis = -6./5;
232         assert(std::abs((mean - x_mean) / x_mean) < 0.01);
233         assert(std::abs((var - x_var) / x_var) < 0.01);
234         assert(std::abs(skew - x_skew) < 0.02);
235         assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.01);
236     }
237     {
238         typedef std::uniform_real_distribution<> D;
239         typedef std::ranlux48_base G;
240         G g;
241         D d;
242         const int N = 100000;
243         std::vector<D::result_type> u;
244         for (int i = 0; i < N; ++i)
245         {
246             D::result_type v = d(g);
247             assert(d.a() <= v && v < d.b());
248             u.push_back(v);
249         }
250         D::result_type mean = std::accumulate(u.begin(), u.end(),
251                                               D::result_type(0)) / u.size();
252         D::result_type var = 0;
253         D::result_type skew = 0;
254         D::result_type kurtosis = 0;
255         for (std::size_t i = 0; i < u.size(); ++i)
256         {
257             D::result_type dbl = (u[i] - mean);
258             D::result_type d2 = sqr(dbl);
259             var += d2;
260             skew += dbl * d2;
261             kurtosis += d2 * d2;
262         }
263         var /= u.size();
264         D::result_type dev = std::sqrt(var);
265         skew /= u.size() * dev * var;
266         kurtosis /= u.size() * var * var;
267         kurtosis -= 3;
268         D::result_type x_mean = (d.a() + d.b()) / 2;
269         D::result_type x_var = sqr(d.b() - d.a()) / 12;
270         D::result_type x_skew = 0;
271         D::result_type x_kurtosis = -6./5;
272         assert(std::abs((mean - x_mean) / x_mean) < 0.01);
273         assert(std::abs((var - x_var) / x_var) < 0.01);
274         assert(std::abs(skew - x_skew) < 0.01);
275         assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.01);
276     }
277     {
278         typedef std::uniform_real_distribution<> D;
279         typedef std::ranlux24 G;
280         G g;
281         D d;
282         const int N = 100000;
283         std::vector<D::result_type> u;
284         for (int i = 0; i < N; ++i)
285         {
286             D::result_type v = d(g);
287             assert(d.a() <= v && v < d.b());
288             u.push_back(v);
289         }
290         D::result_type mean = std::accumulate(u.begin(), u.end(),
291                                               D::result_type(0)) / u.size();
292         D::result_type var = 0;
293         D::result_type skew = 0;
294         D::result_type kurtosis = 0;
295         for (std::size_t i = 0; i < u.size(); ++i)
296         {
297             D::result_type dbl = (u[i] - mean);
298             D::result_type d2 = sqr(dbl);
299             var += d2;
300             skew += dbl * d2;
301             kurtosis += d2 * d2;
302         }
303         var /= u.size();
304         D::result_type dev = std::sqrt(var);
305         skew /= u.size() * dev * var;
306         kurtosis /= u.size() * var * var;
307         kurtosis -= 3;
308         D::result_type x_mean = (d.a() + d.b()) / 2;
309         D::result_type x_var = sqr(d.b() - d.a()) / 12;
310         D::result_type x_skew = 0;
311         D::result_type x_kurtosis = -6./5;
312         assert(std::abs((mean - x_mean) / x_mean) < 0.01);
313         assert(std::abs((var - x_var) / x_var) < 0.01);
314         assert(std::abs(skew - x_skew) < 0.01);
315         assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.01);
316     }
317     {
318         typedef std::uniform_real_distribution<> D;
319         typedef std::ranlux48 G;
320         G g;
321         D d;
322         const int N = 100000;
323         std::vector<D::result_type> u;
324         for (int i = 0; i < N; ++i)
325         {
326             D::result_type v = d(g);
327             assert(d.a() <= v && v < d.b());
328             u.push_back(v);
329         }
330         D::result_type mean = std::accumulate(u.begin(), u.end(),
331                                               D::result_type(0)) / u.size();
332         D::result_type var = 0;
333         D::result_type skew = 0;
334         D::result_type kurtosis = 0;
335         for (std::size_t i = 0; i < u.size(); ++i)
336         {
337             D::result_type dbl = (u[i] - mean);
338             D::result_type d2 = sqr(dbl);
339             var += d2;
340             skew += dbl * d2;
341             kurtosis += d2 * d2;
342         }
343         var /= u.size();
344         D::result_type dev = std::sqrt(var);
345         skew /= u.size() * dev * var;
346         kurtosis /= u.size() * var * var;
347         kurtosis -= 3;
348         D::result_type x_mean = (d.a() + d.b()) / 2;
349         D::result_type x_var = sqr(d.b() - d.a()) / 12;
350         D::result_type x_skew = 0;
351         D::result_type x_kurtosis = -6./5;
352         assert(std::abs((mean - x_mean) / x_mean) < 0.01);
353         assert(std::abs((var - x_var) / x_var) < 0.01);
354         assert(std::abs(skew - x_skew) < 0.01);
355         assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.01);
356     }
357     {
358         typedef std::uniform_real_distribution<> D;
359         typedef std::knuth_b G;
360         G g;
361         D d;
362         const int N = 100000;
363         std::vector<D::result_type> u;
364         for (int i = 0; i < N; ++i)
365         {
366             D::result_type v = d(g);
367             assert(d.a() <= v && v < d.b());
368             u.push_back(v);
369         }
370         D::result_type mean = std::accumulate(u.begin(), u.end(),
371                                               D::result_type(0)) / u.size();
372         D::result_type var = 0;
373         D::result_type skew = 0;
374         D::result_type kurtosis = 0;
375         for (std::size_t i = 0; i < u.size(); ++i)
376         {
377             D::result_type dbl = (u[i] - mean);
378             D::result_type d2 = sqr(dbl);
379             var += d2;
380             skew += dbl * d2;
381             kurtosis += d2 * d2;
382         }
383         var /= u.size();
384         D::result_type dev = std::sqrt(var);
385         skew /= u.size() * dev * var;
386         kurtosis /= u.size() * var * var;
387         kurtosis -= 3;
388         D::result_type x_mean = (d.a() + d.b()) / 2;
389         D::result_type x_var = sqr(d.b() - d.a()) / 12;
390         D::result_type x_skew = 0;
391         D::result_type x_kurtosis = -6./5;
392         assert(std::abs((mean - x_mean) / x_mean) < 0.01);
393         assert(std::abs((var - x_var) / x_var) < 0.01);
394         assert(std::abs(skew - x_skew) < 0.01);
395         assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.01);
396     }
397     {
398         typedef std::uniform_real_distribution<> D;
399         typedef std::minstd_rand G;
400         G g;
401         D d(-1, 1);
402         const int N = 100000;
403         std::vector<D::result_type> u;
404         for (int i = 0; i < N; ++i)
405         {
406             D::result_type v = d(g);
407             assert(d.a() <= v && v < d.b());
408             u.push_back(v);
409         }
410         D::result_type mean = std::accumulate(u.begin(), u.end(),
411                                               D::result_type(0)) / u.size();
412         D::result_type var = 0;
413         D::result_type skew = 0;
414         D::result_type kurtosis = 0;
415         for (std::size_t i = 0; i < u.size(); ++i)
416         {
417             D::result_type dbl = (u[i] - mean);
418             D::result_type d2 = sqr(dbl);
419             var += d2;
420             skew += dbl * d2;
421             kurtosis += d2 * d2;
422         }
423         var /= u.size();
424         D::result_type dev = std::sqrt(var);
425         skew /= u.size() * dev * var;
426         kurtosis /= u.size() * var * var;
427         kurtosis -= 3;
428         D::result_type x_mean = (d.a() + d.b()) / 2;
429         D::result_type x_var = sqr(d.b() - d.a()) / 12;
430         D::result_type x_skew = 0;
431         D::result_type x_kurtosis = -6./5;
432         assert(std::abs(mean - x_mean) < 0.01);
433         assert(std::abs((var - x_var) / x_var) < 0.01);
434         assert(std::abs(skew - x_skew) < 0.01);
435         assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.01);
436     }
437     {
438         typedef std::uniform_real_distribution<> D;
439         typedef std::minstd_rand G;
440         G g;
441         D d(5.5, 25);
442         const int N = 100000;
443         std::vector<D::result_type> u;
444         for (int i = 0; i < N; ++i)
445         {
446             D::result_type v = d(g);
447             assert(d.a() <= v && v < d.b());
448             u.push_back(v);
449         }
450         D::result_type mean = std::accumulate(u.begin(), u.end(),
451                                               D::result_type(0)) / u.size();
452         D::result_type var = 0;
453         D::result_type skew = 0;
454         D::result_type kurtosis = 0;
455         for (std::size_t i = 0; i < u.size(); ++i)
456         {
457             D::result_type dbl = (u[i] - mean);
458             D::result_type d2 = sqr(dbl);
459             var += d2;
460             skew += dbl * d2;
461             kurtosis += d2 * d2;
462         }
463         var /= u.size();
464         D::result_type dev = std::sqrt(var);
465         skew /= u.size() * dev * var;
466         kurtosis /= u.size() * var * var;
467         kurtosis -= 3;
468         D::result_type x_mean = (d.a() + d.b()) / 2;
469         D::result_type x_var = sqr(d.b() - d.a()) / 12;
470         D::result_type x_skew = 0;
471         D::result_type x_kurtosis = -6./5;
472         assert(std::abs((mean - x_mean) / x_mean) < 0.01);
473         assert(std::abs((var - x_var) / x_var) < 0.01);
474         assert(std::abs(skew - x_skew) < 0.01);
475         assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.01);
476     }
477 
478   return 0;
479 }
480