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