1 //===----------------------------------------------------------------------===//
2 //
3 //                     The LLVM Compiler Infrastructure
4 //
5 // This file is dual licensed under the MIT and the University of Illinois Open
6 // Source Licenses. See LICENSE.TXT for details.
7 //
8 //===----------------------------------------------------------------------===//
9 //
10 // REQUIRES: long_tests
11 
12 // <random>
13 
14 // template<class IntType = int>
15 // class geometric_distribution
16 
17 // template<class _URNG> result_type operator()(_URNG& g);
18 
19 #include <random>
20 #include <numeric>
21 #include <vector>
22 #include <cassert>
23 
24 template <class T>
25 inline
26 T
sqr(T x)27 sqr(T x)
28 {
29     return x * x;
30 }
31 
main()32 int main()
33 {
34     {
35         typedef std::geometric_distribution<> D;
36         typedef std::mt19937 G;
37         G g;
38         D d(.03125);
39         const int N = 1000000;
40         std::vector<D::result_type> u;
41         for (int i = 0; i < N; ++i)
42         {
43             D::result_type v = d(g);
44             assert(d.min() <= v && v <= d.max());
45             u.push_back(v);
46         }
47         double mean = std::accumulate(u.begin(), u.end(),
48                                               double(0)) / u.size();
49         double var = 0;
50         double skew = 0;
51         double kurtosis = 0;
52         for (int i = 0; i < u.size(); ++i)
53         {
54             double d = (u[i] - mean);
55             double d2 = sqr(d);
56             var += d2;
57             skew += d * d2;
58             kurtosis += d2 * d2;
59         }
60         var /= u.size();
61         double dev = std::sqrt(var);
62         skew /= u.size() * dev * var;
63         kurtosis /= u.size() * var * var;
64         kurtosis -= 3;
65         double x_mean = (1 - d.p()) / d.p();
66         double x_var = x_mean / d.p();
67         double x_skew = (2 - d.p()) / std::sqrt((1 - d.p()));
68         double x_kurtosis = 6 + sqr(d.p()) / (1 - d.p());
69         assert(std::abs((mean - x_mean) / x_mean) < 0.01);
70         assert(std::abs((var - x_var) / x_var) < 0.01);
71         assert(std::abs((skew - x_skew) / x_skew) < 0.01);
72         assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.01);
73     }
74     {
75         typedef std::geometric_distribution<> D;
76         typedef std::mt19937 G;
77         G g;
78         D d(0.05);
79         const int N = 1000000;
80         std::vector<D::result_type> u;
81         for (int i = 0; i < N; ++i)
82         {
83             D::result_type v = d(g);
84             assert(d.min() <= v && v <= d.max());
85             u.push_back(v);
86         }
87         double mean = std::accumulate(u.begin(), u.end(),
88                                               double(0)) / u.size();
89         double var = 0;
90         double skew = 0;
91         double kurtosis = 0;
92         for (int i = 0; i < u.size(); ++i)
93         {
94             double d = (u[i] - mean);
95             double d2 = sqr(d);
96             var += d2;
97             skew += d * d2;
98             kurtosis += d2 * d2;
99         }
100         var /= u.size();
101         double dev = std::sqrt(var);
102         skew /= u.size() * dev * var;
103         kurtosis /= u.size() * var * var;
104         kurtosis -= 3;
105         double x_mean = (1 - d.p()) / d.p();
106         double x_var = x_mean / d.p();
107         double x_skew = (2 - d.p()) / std::sqrt((1 - d.p()));
108         double x_kurtosis = 6 + sqr(d.p()) / (1 - d.p());
109         assert(std::abs((mean - x_mean) / x_mean) < 0.01);
110         assert(std::abs((var - x_var) / x_var) < 0.01);
111         assert(std::abs((skew - x_skew) / x_skew) < 0.01);
112         assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.03);
113     }
114     {
115         typedef std::geometric_distribution<> D;
116         typedef std::minstd_rand G;
117         G g;
118         D d(.25);
119         const int N = 1000000;
120         std::vector<D::result_type> u;
121         for (int i = 0; i < N; ++i)
122         {
123             D::result_type v = d(g);
124             assert(d.min() <= v && v <= d.max());
125             u.push_back(v);
126         }
127         double mean = std::accumulate(u.begin(), u.end(),
128                                               double(0)) / u.size();
129         double var = 0;
130         double skew = 0;
131         double kurtosis = 0;
132         for (int i = 0; i < u.size(); ++i)
133         {
134             double d = (u[i] - mean);
135             double d2 = sqr(d);
136             var += d2;
137             skew += d * d2;
138             kurtosis += d2 * d2;
139         }
140         var /= u.size();
141         double dev = std::sqrt(var);
142         skew /= u.size() * dev * var;
143         kurtosis /= u.size() * var * var;
144         kurtosis -= 3;
145         double x_mean = (1 - d.p()) / d.p();
146         double x_var = x_mean / d.p();
147         double x_skew = (2 - d.p()) / std::sqrt((1 - d.p()));
148         double x_kurtosis = 6 + sqr(d.p()) / (1 - d.p());
149         assert(std::abs((mean - x_mean) / x_mean) < 0.01);
150         assert(std::abs((var - x_var) / x_var) < 0.01);
151         assert(std::abs((skew - x_skew) / x_skew) < 0.01);
152         assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.02);
153     }
154     {
155         typedef std::geometric_distribution<> D;
156         typedef std::mt19937 G;
157         G g;
158         D d(0.5);
159         const int N = 1000000;
160         std::vector<D::result_type> u;
161         for (int i = 0; i < N; ++i)
162         {
163             D::result_type v = d(g);
164             assert(d.min() <= v && v <= d.max());
165             u.push_back(v);
166         }
167         double mean = std::accumulate(u.begin(), u.end(),
168                                               double(0)) / u.size();
169         double var = 0;
170         double skew = 0;
171         double kurtosis = 0;
172         for (int i = 0; i < u.size(); ++i)
173         {
174             double d = (u[i] - mean);
175             double d2 = sqr(d);
176             var += d2;
177             skew += d * d2;
178             kurtosis += d2 * d2;
179         }
180         var /= u.size();
181         double dev = std::sqrt(var);
182         skew /= u.size() * dev * var;
183         kurtosis /= u.size() * var * var;
184         kurtosis -= 3;
185         double x_mean = (1 - d.p()) / d.p();
186         double x_var = x_mean / d.p();
187         double x_skew = (2 - d.p()) / std::sqrt((1 - d.p()));
188         double x_kurtosis = 6 + sqr(d.p()) / (1 - d.p());
189         assert(std::abs((mean - x_mean) / x_mean) < 0.01);
190         assert(std::abs((var - x_var) / x_var) < 0.01);
191         assert(std::abs((skew - x_skew) / x_skew) < 0.01);
192         assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.02);
193     }
194     {
195         typedef std::geometric_distribution<> D;
196         typedef std::mt19937 G;
197         G g;
198         D d(0.75);
199         const int N = 1000000;
200         std::vector<D::result_type> u;
201         for (int i = 0; i < N; ++i)
202         {
203             D::result_type v = d(g);
204             assert(d.min() <= v && v <= d.max());
205             u.push_back(v);
206         }
207         double mean = std::accumulate(u.begin(), u.end(),
208                                               double(0)) / u.size();
209         double var = 0;
210         double skew = 0;
211         double kurtosis = 0;
212         for (int i = 0; i < u.size(); ++i)
213         {
214             double d = (u[i] - mean);
215             double d2 = sqr(d);
216             var += d2;
217             skew += d * d2;
218             kurtosis += d2 * d2;
219         }
220         var /= u.size();
221         double dev = std::sqrt(var);
222         skew /= u.size() * dev * var;
223         kurtosis /= u.size() * var * var;
224         kurtosis -= 3;
225         double x_mean = (1 - d.p()) / d.p();
226         double x_var = x_mean / d.p();
227         double x_skew = (2 - d.p()) / std::sqrt((1 - d.p()));
228         double x_kurtosis = 6 + sqr(d.p()) / (1 - d.p());
229         assert(std::abs((mean - x_mean) / x_mean) < 0.01);
230         assert(std::abs((var - x_var) / x_var) < 0.01);
231         assert(std::abs((skew - x_skew) / x_skew) < 0.01);
232         assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.02);
233     }
234     {
235         typedef std::geometric_distribution<> D;
236         typedef std::mt19937 G;
237         G g;
238         D d(0.96875);
239         const int N = 1000000;
240         std::vector<D::result_type> u;
241         for (int i = 0; i < N; ++i)
242         {
243             D::result_type v = d(g);
244             assert(d.min() <= v && v <= d.max());
245             u.push_back(v);
246         }
247         double mean = std::accumulate(u.begin(), u.end(),
248                                               double(0)) / u.size();
249         double var = 0;
250         double skew = 0;
251         double kurtosis = 0;
252         for (int i = 0; i < u.size(); ++i)
253         {
254             double d = (u[i] - mean);
255             double d2 = sqr(d);
256             var += d2;
257             skew += d * d2;
258             kurtosis += d2 * d2;
259         }
260         var /= u.size();
261         double dev = std::sqrt(var);
262         skew /= u.size() * dev * var;
263         kurtosis /= u.size() * var * var;
264         kurtosis -= 3;
265         double x_mean = (1 - d.p()) / d.p();
266         double x_var = x_mean / d.p();
267         double x_skew = (2 - d.p()) / std::sqrt((1 - d.p()));
268         double x_kurtosis = 6 + sqr(d.p()) / (1 - d.p());
269         assert(std::abs((mean - x_mean) / x_mean) < 0.01);
270         assert(std::abs((var - x_var) / x_var) < 0.01);
271         assert(std::abs((skew - x_skew) / x_skew) < 0.01);
272         assert(std::abs((kurtosis - x_kurtosis) / x_kurtosis) < 0.02);
273     }
274 }
275