libstdc++
bits/random.tcc
Go to the documentation of this file.
1// random number generation (out of line) -*- C++ -*-
2
3// Copyright (C) 2009-2025 Free Software Foundation, Inc.
4//
5// This file is part of the GNU ISO C++ Library. This library is free
6// software; you can redistribute it and/or modify it under the
7// terms of the GNU General Public License as published by the
8// Free Software Foundation; either version 3, or (at your option)
9// any later version.
10
11// This library is distributed in the hope that it will be useful,
12// but WITHOUT ANY WARRANTY; without even the implied warranty of
13// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14// GNU General Public License for more details.
15
16// Under Section 7 of GPL version 3, you are granted additional
17// permissions described in the GCC Runtime Library Exception, version
18// 3.1, as published by the Free Software Foundation.
19
20// You should have received a copy of the GNU General Public License and
21// a copy of the GCC Runtime Library Exception along with this program;
22// see the files COPYING3 and COPYING.RUNTIME respectively. If not, see
23// <http://www.gnu.org/licenses/>.
24
25/** @file bits/random.tcc
26 * This is an internal header file, included by other library headers.
27 * Do not attempt to use it directly. @headername{random}
28 */
29
30#ifndef _RANDOM_TCC
31#define _RANDOM_TCC 1
32
33#include <numeric> // std::accumulate and std::partial_sum
34
35namespace std _GLIBCXX_VISIBILITY(default)
36{
37_GLIBCXX_BEGIN_NAMESPACE_VERSION
38
39 /// @cond undocumented
40 // (Further) implementation-space details.
41 namespace __detail
42 {
43 // General case for x = (ax + c) mod m -- use Schrage's algorithm
44 // to avoid integer overflow.
45 //
46 // Preconditions: a > 0, m > 0.
47 //
48 // Note: only works correctly for __m % __a < __m / __a.
49 template<typename _Tp, _Tp __m, _Tp __a, _Tp __c>
50 _Tp
51 _Mod<_Tp, __m, __a, __c, false, true>::
52 __calc(_Tp __x)
53 {
54 if (__a == 1)
55 __x %= __m;
56 else
57 {
58 static const _Tp __q = __m / __a;
59 static const _Tp __r = __m % __a;
60
61 _Tp __t1 = __a * (__x % __q);
62 _Tp __t2 = __r * (__x / __q);
63 if (__t1 >= __t2)
64 __x = __t1 - __t2;
65 else
66 __x = __m - __t2 + __t1;
67 }
68
69 if (__c != 0)
70 {
71 const _Tp __d = __m - __x;
72 if (__d > __c)
73 __x += __c;
74 else
75 __x = __c - __d;
76 }
77 return __x;
78 }
79
80 template<typename _InputIterator, typename _OutputIterator,
81 typename _Tp>
82 _OutputIterator
83 __normalize(_InputIterator __first, _InputIterator __last,
84 _OutputIterator __result, const _Tp& __factor)
85 {
86 for (; __first != __last; ++__first, (void) ++__result)
87 *__result = *__first / __factor;
88 return __result;
89 }
90
91 } // namespace __detail
92 /// @endcond
93
94#if ! __cpp_inline_variables
95 template<typename _UIntType, _UIntType __a, _UIntType __c, _UIntType __m>
96 constexpr _UIntType
98
99 template<typename _UIntType, _UIntType __a, _UIntType __c, _UIntType __m>
100 constexpr _UIntType
102
103 template<typename _UIntType, _UIntType __a, _UIntType __c, _UIntType __m>
104 constexpr _UIntType
106
107 template<typename _UIntType, _UIntType __a, _UIntType __c, _UIntType __m>
108 constexpr _UIntType
109 linear_congruential_engine<_UIntType, __a, __c, __m>::default_seed;
110#endif
111
112 /**
113 * Seeds the LCR with integral value @p __s, adjusted so that the
114 * ring identity is never a member of the convergence set.
115 */
116 template<typename _UIntType, _UIntType __a, _UIntType __c, _UIntType __m>
117 void
120 {
121 if ((__detail::__mod<_UIntType, __m>(__c) == 0)
122 && (__detail::__mod<_UIntType, __m>(__s) == 0))
123 _M_x = 1;
124 else
125 _M_x = __detail::__mod<_UIntType, __m>(__s);
126 }
127
128 /**
129 * Seeds the LCR engine with a value generated by @p __q.
130 */
131 template<typename _UIntType, _UIntType __a, _UIntType __c, _UIntType __m>
132 template<typename _Sseq>
133 auto
135 seed(_Sseq& __q)
136 -> _If_seed_seq<_Sseq>
137 {
138 const _UIntType __k0 = __m == 0 ? std::numeric_limits<_UIntType>::digits
139 : std::__lg(__m);
140 const _UIntType __k = (__k0 + 31) / 32;
141 uint_least32_t __arr[__k + 3];
142 __q.generate(__arr + 0, __arr + __k + 3);
143 _UIntType __factor = 1u;
144 _UIntType __sum = 0u;
145 for (size_t __j = 0; __j < __k; ++__j)
146 {
147 __sum += __arr[__j + 3] * __factor;
148 __factor *= __detail::_Shift<_UIntType, 32>::__value;
149 }
150 seed(__sum);
151 }
152
153 template<typename _UIntType, _UIntType __a, _UIntType __c, _UIntType __m,
154 typename _CharT, typename _Traits>
157 const linear_congruential_engine<_UIntType,
158 __a, __c, __m>& __lcr)
159 {
160 using __ios_base = typename basic_ostream<_CharT, _Traits>::ios_base;
161
162 const typename __ios_base::fmtflags __flags = __os.flags();
163 const _CharT __fill = __os.fill();
164 __os.flags(__ios_base::dec | __ios_base::fixed | __ios_base::left);
165 __os.fill(__os.widen(' '));
166
167 __os << __lcr._M_x;
168
169 __os.flags(__flags);
170 __os.fill(__fill);
171 return __os;
172 }
173
174 template<typename _UIntType, _UIntType __a, _UIntType __c, _UIntType __m,
175 typename _CharT, typename _Traits>
178 linear_congruential_engine<_UIntType, __a, __c, __m>& __lcr)
179 {
180 using __ios_base = typename basic_istream<_CharT, _Traits>::ios_base;
181
182 const typename __ios_base::fmtflags __flags = __is.flags();
183 __is.flags(__ios_base::dec);
184
185 __is >> __lcr._M_x;
186
187 __is.flags(__flags);
188 return __is;
189 }
190
191#if ! __cpp_inline_variables
192 template<typename _UIntType,
193 size_t __w, size_t __n, size_t __m, size_t __r,
194 _UIntType __a, size_t __u, _UIntType __d, size_t __s,
195 _UIntType __b, size_t __t, _UIntType __c, size_t __l,
196 _UIntType __f>
197 constexpr size_t
198 mersenne_twister_engine<_UIntType, __w, __n, __m, __r, __a, __u, __d,
199 __s, __b, __t, __c, __l, __f>::word_size;
200
201 template<typename _UIntType,
202 size_t __w, size_t __n, size_t __m, size_t __r,
203 _UIntType __a, size_t __u, _UIntType __d, size_t __s,
204 _UIntType __b, size_t __t, _UIntType __c, size_t __l,
205 _UIntType __f>
206 constexpr size_t
207 mersenne_twister_engine<_UIntType, __w, __n, __m, __r, __a, __u, __d,
208 __s, __b, __t, __c, __l, __f>::state_size;
209
210 template<typename _UIntType,
211 size_t __w, size_t __n, size_t __m, size_t __r,
212 _UIntType __a, size_t __u, _UIntType __d, size_t __s,
213 _UIntType __b, size_t __t, _UIntType __c, size_t __l,
214 _UIntType __f>
215 constexpr size_t
216 mersenne_twister_engine<_UIntType, __w, __n, __m, __r, __a, __u, __d,
217 __s, __b, __t, __c, __l, __f>::shift_size;
218
219 template<typename _UIntType,
220 size_t __w, size_t __n, size_t __m, size_t __r,
221 _UIntType __a, size_t __u, _UIntType __d, size_t __s,
222 _UIntType __b, size_t __t, _UIntType __c, size_t __l,
223 _UIntType __f>
224 constexpr size_t
225 mersenne_twister_engine<_UIntType, __w, __n, __m, __r, __a, __u, __d,
226 __s, __b, __t, __c, __l, __f>::mask_bits;
227
228 template<typename _UIntType,
229 size_t __w, size_t __n, size_t __m, size_t __r,
230 _UIntType __a, size_t __u, _UIntType __d, size_t __s,
231 _UIntType __b, size_t __t, _UIntType __c, size_t __l,
232 _UIntType __f>
233 constexpr _UIntType
234 mersenne_twister_engine<_UIntType, __w, __n, __m, __r, __a, __u, __d,
235 __s, __b, __t, __c, __l, __f>::xor_mask;
236
237 template<typename _UIntType,
238 size_t __w, size_t __n, size_t __m, size_t __r,
239 _UIntType __a, size_t __u, _UIntType __d, size_t __s,
240 _UIntType __b, size_t __t, _UIntType __c, size_t __l,
241 _UIntType __f>
242 constexpr size_t
243 mersenne_twister_engine<_UIntType, __w, __n, __m, __r, __a, __u, __d,
244 __s, __b, __t, __c, __l, __f>::tempering_u;
245
246 template<typename _UIntType,
247 size_t __w, size_t __n, size_t __m, size_t __r,
248 _UIntType __a, size_t __u, _UIntType __d, size_t __s,
249 _UIntType __b, size_t __t, _UIntType __c, size_t __l,
250 _UIntType __f>
251 constexpr _UIntType
252 mersenne_twister_engine<_UIntType, __w, __n, __m, __r, __a, __u, __d,
253 __s, __b, __t, __c, __l, __f>::tempering_d;
254
255 template<typename _UIntType,
256 size_t __w, size_t __n, size_t __m, size_t __r,
257 _UIntType __a, size_t __u, _UIntType __d, size_t __s,
258 _UIntType __b, size_t __t, _UIntType __c, size_t __l,
259 _UIntType __f>
260 constexpr size_t
261 mersenne_twister_engine<_UIntType, __w, __n, __m, __r, __a, __u, __d,
262 __s, __b, __t, __c, __l, __f>::tempering_s;
263
264 template<typename _UIntType,
265 size_t __w, size_t __n, size_t __m, size_t __r,
266 _UIntType __a, size_t __u, _UIntType __d, size_t __s,
267 _UIntType __b, size_t __t, _UIntType __c, size_t __l,
268 _UIntType __f>
269 constexpr _UIntType
270 mersenne_twister_engine<_UIntType, __w, __n, __m, __r, __a, __u, __d,
271 __s, __b, __t, __c, __l, __f>::tempering_b;
272
273 template<typename _UIntType,
274 size_t __w, size_t __n, size_t __m, size_t __r,
275 _UIntType __a, size_t __u, _UIntType __d, size_t __s,
276 _UIntType __b, size_t __t, _UIntType __c, size_t __l,
277 _UIntType __f>
278 constexpr size_t
279 mersenne_twister_engine<_UIntType, __w, __n, __m, __r, __a, __u, __d,
280 __s, __b, __t, __c, __l, __f>::tempering_t;
281
282 template<typename _UIntType,
283 size_t __w, size_t __n, size_t __m, size_t __r,
284 _UIntType __a, size_t __u, _UIntType __d, size_t __s,
285 _UIntType __b, size_t __t, _UIntType __c, size_t __l,
286 _UIntType __f>
287 constexpr _UIntType
288 mersenne_twister_engine<_UIntType, __w, __n, __m, __r, __a, __u, __d,
289 __s, __b, __t, __c, __l, __f>::tempering_c;
290
291 template<typename _UIntType,
292 size_t __w, size_t __n, size_t __m, size_t __r,
293 _UIntType __a, size_t __u, _UIntType __d, size_t __s,
294 _UIntType __b, size_t __t, _UIntType __c, size_t __l,
295 _UIntType __f>
296 constexpr size_t
297 mersenne_twister_engine<_UIntType, __w, __n, __m, __r, __a, __u, __d,
298 __s, __b, __t, __c, __l, __f>::tempering_l;
299
300 template<typename _UIntType,
301 size_t __w, size_t __n, size_t __m, size_t __r,
302 _UIntType __a, size_t __u, _UIntType __d, size_t __s,
303 _UIntType __b, size_t __t, _UIntType __c, size_t __l,
304 _UIntType __f>
305 constexpr _UIntType
306 mersenne_twister_engine<_UIntType, __w, __n, __m, __r, __a, __u, __d,
307 __s, __b, __t, __c, __l, __f>::
308 initialization_multiplier;
309
310 template<typename _UIntType,
311 size_t __w, size_t __n, size_t __m, size_t __r,
312 _UIntType __a, size_t __u, _UIntType __d, size_t __s,
313 _UIntType __b, size_t __t, _UIntType __c, size_t __l,
314 _UIntType __f>
315 constexpr _UIntType
316 mersenne_twister_engine<_UIntType, __w, __n, __m, __r, __a, __u, __d,
317 __s, __b, __t, __c, __l, __f>::default_seed;
318#endif
319
320 template<typename _UIntType,
321 size_t __w, size_t __n, size_t __m, size_t __r,
322 _UIntType __a, size_t __u, _UIntType __d, size_t __s,
323 _UIntType __b, size_t __t, _UIntType __c, size_t __l,
324 _UIntType __f>
325 void
326 mersenne_twister_engine<_UIntType, __w, __n, __m, __r, __a, __u, __d,
327 __s, __b, __t, __c, __l, __f>::
328 seed(result_type __sd)
329 {
330 _M_x[0] = __detail::__mod<_UIntType,
331 __detail::_Shift<_UIntType, __w>::__value>(__sd);
332
333 for (size_t __i = 1; __i < state_size; ++__i)
334 {
335 _UIntType __x = _M_x[__i - 1];
336 __x ^= __x >> (__w - 2);
337 __x *= __f;
338 __x += __detail::__mod<_UIntType, __n>(__i);
339 _M_x[__i] = __detail::__mod<_UIntType,
340 __detail::_Shift<_UIntType, __w>::__value>(__x);
341 }
342 _M_p = state_size;
343 }
344
345 template<typename _UIntType,
346 size_t __w, size_t __n, size_t __m, size_t __r,
347 _UIntType __a, size_t __u, _UIntType __d, size_t __s,
348 _UIntType __b, size_t __t, _UIntType __c, size_t __l,
349 _UIntType __f>
350 template<typename _Sseq>
351 auto
352 mersenne_twister_engine<_UIntType, __w, __n, __m, __r, __a, __u, __d,
353 __s, __b, __t, __c, __l, __f>::
354 seed(_Sseq& __q)
355 -> _If_seed_seq<_Sseq>
356 {
357 const _UIntType __upper_mask = (~_UIntType()) << __r;
358 const size_t __k = (__w + 31) / 32;
359 uint_least32_t __arr[__n * __k];
360 __q.generate(__arr + 0, __arr + __n * __k);
361
362 bool __zero = true;
363 for (size_t __i = 0; __i < state_size; ++__i)
364 {
365 _UIntType __factor = 1u;
366 _UIntType __sum = 0u;
367 for (size_t __j = 0; __j < __k; ++__j)
368 {
369 __sum += __arr[__k * __i + __j] * __factor;
370 __factor *= __detail::_Shift<_UIntType, 32>::__value;
371 }
372 _M_x[__i] = __detail::__mod<_UIntType,
373 __detail::_Shift<_UIntType, __w>::__value>(__sum);
374
375 if (__zero)
376 {
377 if (__i == 0)
378 {
379 if ((_M_x[0] & __upper_mask) != 0u)
380 __zero = false;
381 }
382 else if (_M_x[__i] != 0u)
383 __zero = false;
384 }
385 }
386 if (__zero)
387 _M_x[0] = __detail::_Shift<_UIntType, __w - 1>::__value;
388 _M_p = state_size;
389 }
390
391 template<typename _UIntType, size_t __w,
392 size_t __n, size_t __m, size_t __r,
393 _UIntType __a, size_t __u, _UIntType __d, size_t __s,
394 _UIntType __b, size_t __t, _UIntType __c, size_t __l,
395 _UIntType __f>
396 void
397 mersenne_twister_engine<_UIntType, __w, __n, __m, __r, __a, __u, __d,
398 __s, __b, __t, __c, __l, __f>::
399 _M_gen_rand(void)
400 {
401 const _UIntType __upper_mask = (~_UIntType()) << __r;
402 const _UIntType __lower_mask = ~__upper_mask;
403
404 for (size_t __k = 0; __k < (__n - __m); ++__k)
405 {
406 _UIntType __y = ((_M_x[__k] & __upper_mask)
407 | (_M_x[__k + 1] & __lower_mask));
408 _M_x[__k] = (_M_x[__k + __m] ^ (__y >> 1)
409 ^ ((__y & 0x01) ? __a : 0));
410 }
411
412 for (size_t __k = (__n - __m); __k < (__n - 1); ++__k)
413 {
414 _UIntType __y = ((_M_x[__k] & __upper_mask)
415 | (_M_x[__k + 1] & __lower_mask));
416 _M_x[__k] = (_M_x[__k + (__m - __n)] ^ (__y >> 1)
417 ^ ((__y & 0x01) ? __a : 0));
418 }
419
420 _UIntType __y = ((_M_x[__n - 1] & __upper_mask)
421 | (_M_x[0] & __lower_mask));
422 _M_x[__n - 1] = (_M_x[__m - 1] ^ (__y >> 1)
423 ^ ((__y & 0x01) ? __a : 0));
424 _M_p = 0;
425 }
427 template<typename _UIntType, size_t __w,
428 size_t __n, size_t __m, size_t __r,
429 _UIntType __a, size_t __u, _UIntType __d, size_t __s,
430 _UIntType __b, size_t __t, _UIntType __c, size_t __l,
431 _UIntType __f>
432 void
433 mersenne_twister_engine<_UIntType, __w, __n, __m, __r, __a, __u, __d,
434 __s, __b, __t, __c, __l, __f>::
435 discard(unsigned long long __z)
436 {
437 while (__z > state_size - _M_p)
438 {
439 __z -= state_size - _M_p;
440 _M_gen_rand();
441 }
442 _M_p += __z;
443 }
444
445 template<typename _UIntType, size_t __w,
446 size_t __n, size_t __m, size_t __r,
447 _UIntType __a, size_t __u, _UIntType __d, size_t __s,
448 _UIntType __b, size_t __t, _UIntType __c, size_t __l,
449 _UIntType __f>
450 typename
451 mersenne_twister_engine<_UIntType, __w, __n, __m, __r, __a, __u, __d,
452 __s, __b, __t, __c, __l, __f>::result_type
453 mersenne_twister_engine<_UIntType, __w, __n, __m, __r, __a, __u, __d,
454 __s, __b, __t, __c, __l, __f>::
455 operator()()
456 {
457 // Reload the vector - cost is O(n) amortized over n calls.
458 if (_M_p >= state_size)
459 _M_gen_rand();
460
461 // Calculate o(x(i)).
462 result_type __z = _M_x[_M_p++];
463 __z ^= (__z >> __u) & __d;
464 __z ^= (__z << __s) & __b;
465 __z ^= (__z << __t) & __c;
466 __z ^= (__z >> __l);
467
468 return __z;
469 }
470
471 template<typename _UIntType, size_t __w,
472 size_t __n, size_t __m, size_t __r,
473 _UIntType __a, size_t __u, _UIntType __d, size_t __s,
474 _UIntType __b, size_t __t, _UIntType __c, size_t __l,
475 _UIntType __f, typename _CharT, typename _Traits>
478 const mersenne_twister_engine<_UIntType, __w, __n, __m,
479 __r, __a, __u, __d, __s, __b, __t, __c, __l, __f>& __x)
480 {
481 using __ios_base = typename basic_ostream<_CharT, _Traits>::ios_base;
482
483 const typename __ios_base::fmtflags __flags = __os.flags();
484 const _CharT __fill = __os.fill();
485 const _CharT __space = __os.widen(' ');
486 __os.flags(__ios_base::dec | __ios_base::fixed | __ios_base::left);
487 __os.fill(__space);
488
489 for (size_t __i = 0; __i < __n; ++__i)
490 __os << __x._M_x[__i] << __space;
491 __os << __x._M_p;
492
493 __os.flags(__flags);
494 __os.fill(__fill);
495 return __os;
496 }
497
498 template<typename _UIntType, size_t __w,
499 size_t __n, size_t __m, size_t __r,
500 _UIntType __a, size_t __u, _UIntType __d, size_t __s,
501 _UIntType __b, size_t __t, _UIntType __c, size_t __l,
502 _UIntType __f, typename _CharT, typename _Traits>
503 std::basic_istream<_CharT, _Traits>&
504 operator>>(std::basic_istream<_CharT, _Traits>& __is,
505 mersenne_twister_engine<_UIntType, __w, __n, __m,
506 __r, __a, __u, __d, __s, __b, __t, __c, __l, __f>& __x)
507 {
508 using __ios_base = typename basic_istream<_CharT, _Traits>::ios_base;
509
510 const typename __ios_base::fmtflags __flags = __is.flags();
511 __is.flags(__ios_base::dec | __ios_base::skipws);
512
513 for (size_t __i = 0; __i < __n; ++__i)
514 __is >> __x._M_x[__i];
515 __is >> __x._M_p;
516
517 __is.flags(__flags);
518 return __is;
519 }
520
521#if ! __cpp_inline_variables
522 template<typename _UIntType, size_t __w, size_t __s, size_t __r>
523 constexpr size_t
524 subtract_with_carry_engine<_UIntType, __w, __s, __r>::word_size;
525
526 template<typename _UIntType, size_t __w, size_t __s, size_t __r>
527 constexpr size_t
528 subtract_with_carry_engine<_UIntType, __w, __s, __r>::short_lag;
529
530 template<typename _UIntType, size_t __w, size_t __s, size_t __r>
531 constexpr size_t
532 subtract_with_carry_engine<_UIntType, __w, __s, __r>::long_lag;
533
534 template<typename _UIntType, size_t __w, size_t __s, size_t __r>
535 constexpr uint_least32_t
536 subtract_with_carry_engine<_UIntType, __w, __s, __r>::default_seed;
537#endif
538
539 template<typename _UIntType, size_t __w, size_t __s, size_t __r>
540 void
542 seed(result_type __value)
543 {
544 // _GLIBCXX_RESOLVE_LIB_DEFECTS
545 // 3809. Is std::subtract_with_carry_engine<uint16_t> supposed to work?
546 // 4014. LWG 3809 changes behavior of some existing code
548 __lcg(__value == 0u ? default_seed : __value % 2147483563u);
549
550 const size_t __n = (__w + 31) / 32;
551
552 for (size_t __i = 0; __i < long_lag; ++__i)
553 {
554 _UIntType __sum = 0u;
555 _UIntType __factor = 1u;
556 for (size_t __j = 0; __j < __n; ++__j)
557 {
558 __sum += __detail::__mod<uint_least32_t,
559 __detail::_Shift<uint_least32_t, 32>::__value>
560 (__lcg()) * __factor;
561 __factor *= __detail::_Shift<_UIntType, 32>::__value;
562 }
563 _M_x[__i] = __detail::__mod<_UIntType,
564 __detail::_Shift<_UIntType, __w>::__value>(__sum);
565 }
566 _M_carry = (_M_x[long_lag - 1] == 0) ? 1 : 0;
567 _M_p = 0;
568 }
569
570 template<typename _UIntType, size_t __w, size_t __s, size_t __r>
571 template<typename _Sseq>
572 auto
574 seed(_Sseq& __q)
575 -> _If_seed_seq<_Sseq>
576 {
577 const size_t __k = (__w + 31) / 32;
578 uint_least32_t __arr[__r * __k];
579 __q.generate(__arr + 0, __arr + __r * __k);
580
581 for (size_t __i = 0; __i < long_lag; ++__i)
582 {
583 _UIntType __sum = 0u;
584 _UIntType __factor = 1u;
585 for (size_t __j = 0; __j < __k; ++__j)
586 {
587 __sum += __arr[__k * __i + __j] * __factor;
588 __factor *= __detail::_Shift<_UIntType, 32>::__value;
589 }
590 _M_x[__i] = __detail::__mod<_UIntType,
591 __detail::_Shift<_UIntType, __w>::__value>(__sum);
592 }
593 _M_carry = (_M_x[long_lag - 1] == 0) ? 1 : 0;
594 _M_p = 0;
595 }
596
597 template<typename _UIntType, size_t __w, size_t __s, size_t __r>
599 result_type
602 {
603 // Derive short lag index from current index.
604 long __ps = _M_p - short_lag;
605 if (__ps < 0)
606 __ps += long_lag;
607
608 // Calculate new x(i) without overflow or division.
609 // NB: Thanks to the requirements for _UIntType, _M_x[_M_p] + _M_carry
610 // cannot overflow.
611 _UIntType __xi;
612 if (_M_x[__ps] >= _M_x[_M_p] + _M_carry)
613 {
614 __xi = _M_x[__ps] - _M_x[_M_p] - _M_carry;
615 _M_carry = 0;
616 }
617 else
618 {
619 __xi = (__detail::_Shift<_UIntType, __w>::__value
620 - _M_x[_M_p] - _M_carry + _M_x[__ps]);
621 _M_carry = 1;
622 }
623 _M_x[_M_p] = __xi;
624
625 // Adjust current index to loop around in ring buffer.
626 if (++_M_p >= long_lag)
627 _M_p = 0;
628
629 return __xi;
630 }
631
632 template<typename _UIntType, size_t __w, size_t __s, size_t __r,
633 typename _CharT, typename _Traits>
636 const subtract_with_carry_engine<_UIntType,
637 __w, __s, __r>& __x)
638 {
639 using __ios_base = typename basic_ostream<_CharT, _Traits>::ios_base;
640
641 const typename __ios_base::fmtflags __flags = __os.flags();
642 const _CharT __fill = __os.fill();
643 const _CharT __space = __os.widen(' ');
644 __os.flags(__ios_base::dec | __ios_base::fixed | __ios_base::left);
645 __os.fill(__space);
646
647 for (size_t __i = 0; __i < __r; ++__i)
648 __os << __x._M_x[__i] << __space;
649 __os << __x._M_carry << __space << __x._M_p;
650
651 __os.flags(__flags);
652 __os.fill(__fill);
653 return __os;
654 }
655
656 template<typename _UIntType, size_t __w, size_t __s, size_t __r,
657 typename _CharT, typename _Traits>
660 subtract_with_carry_engine<_UIntType, __w, __s, __r>& __x)
661 {
662 using __ios_base = typename basic_istream<_CharT, _Traits>::ios_base;
663
664 const typename __ios_base::fmtflags __flags = __is.flags();
665 __is.flags(__ios_base::dec | __ios_base::skipws);
666
667 for (size_t __i = 0; __i < __r; ++__i)
668 __is >> __x._M_x[__i];
669 __is >> __x._M_carry;
670 __is >> __x._M_p;
671
672 __is.flags(__flags);
673 return __is;
674 }
675
676#if ! __cpp_inline_variables
677 template<typename _RandomNumberEngine, size_t __p, size_t __r>
678 constexpr size_t
679 discard_block_engine<_RandomNumberEngine, __p, __r>::block_size;
680
681 template<typename _RandomNumberEngine, size_t __p, size_t __r>
682 constexpr size_t
683 discard_block_engine<_RandomNumberEngine, __p, __r>::used_block;
684#endif
685
686 template<typename _RandomNumberEngine, size_t __p, size_t __r>
687 typename discard_block_engine<_RandomNumberEngine,
688 __p, __r>::result_type
691 {
692 if (_M_n >= used_block)
693 {
694 _M_b.discard(block_size - _M_n);
695 _M_n = 0;
696 }
697 ++_M_n;
698 return _M_b();
699 }
700
701 template<typename _RandomNumberEngine, size_t __p, size_t __r,
702 typename _CharT, typename _Traits>
705 const discard_block_engine<_RandomNumberEngine,
706 __p, __r>& __x)
707 {
708 using __ios_base = typename basic_ostream<_CharT, _Traits>::ios_base;
709
710 const typename __ios_base::fmtflags __flags = __os.flags();
711 const _CharT __fill = __os.fill();
712 const _CharT __space = __os.widen(' ');
713 __os.flags(__ios_base::dec | __ios_base::fixed | __ios_base::left);
714 __os.fill(__space);
715
716 __os << __x.base() << __space << __x._M_n;
717
718 __os.flags(__flags);
719 __os.fill(__fill);
720 return __os;
721 }
722
723 template<typename _RandomNumberEngine, size_t __p, size_t __r,
724 typename _CharT, typename _Traits>
727 discard_block_engine<_RandomNumberEngine, __p, __r>& __x)
728 {
729 using __ios_base = typename basic_istream<_CharT, _Traits>::ios_base;
730
731 const typename __ios_base::fmtflags __flags = __is.flags();
732 __is.flags(__ios_base::dec | __ios_base::skipws);
733
734 __is >> __x._M_b >> __x._M_n;
735
736 __is.flags(__flags);
737 return __is;
738 }
739
740
741 template<typename _RandomNumberEngine, size_t __w, typename _UIntType>
743 result_type
746 {
747 typedef typename _RandomNumberEngine::result_type _Eresult_type;
748 const _Eresult_type __r
749 = (_M_b.max() - _M_b.min() < std::numeric_limits<_Eresult_type>::max()
750 ? _M_b.max() - _M_b.min() + 1 : 0);
751 const unsigned __edig = std::numeric_limits<_Eresult_type>::digits;
752 const unsigned __m = __r ? std::__lg(__r) : __edig;
753
755 __ctype;
756 const unsigned __cdig = std::numeric_limits<__ctype>::digits;
757
758 unsigned __n, __n0;
759 __ctype __s0, __s1, __y0, __y1;
760
761 for (size_t __i = 0; __i < 2; ++__i)
762 {
763 __n = (__w + __m - 1) / __m + __i;
764 __n0 = __n - __w % __n;
765 const unsigned __w0 = __w / __n; // __w0 <= __m
766
767 __s0 = 0;
768 __s1 = 0;
769 if (__w0 < __cdig)
770 {
771 __s0 = __ctype(1) << __w0;
772 __s1 = __s0 << 1;
773 }
774
775 __y0 = 0;
776 __y1 = 0;
777 if (__r)
778 {
779 __y0 = __s0 * (__r / __s0);
780 if (__s1)
781 __y1 = __s1 * (__r / __s1);
782
783 if (__r - __y0 <= __y0 / __n)
784 break;
785 }
786 else
787 break;
788 }
789
790 result_type __sum = 0;
791 for (size_t __k = 0; __k < __n0; ++__k)
792 {
793 __ctype __u;
794 do
795 __u = _M_b() - _M_b.min();
796 while (__y0 && __u >= __y0);
797 __sum = __s0 * __sum + (__s0 ? __u % __s0 : __u);
798 }
799 for (size_t __k = __n0; __k < __n; ++__k)
800 {
801 __ctype __u;
802 do
803 __u = _M_b() - _M_b.min();
804 while (__y1 && __u >= __y1);
805 __sum = __s1 * __sum + (__s1 ? __u % __s1 : __u);
806 }
807 return __sum;
808 }
809
810#if ! __cpp_inline_variables
811 template<typename _RandomNumberEngine, size_t __k>
812 constexpr size_t
813 shuffle_order_engine<_RandomNumberEngine, __k>::table_size;
814#endif
815
816 namespace __detail
817 {
818 // Determine whether an integer is representable as double.
819 template<typename _Tp>
820 constexpr bool
821 __representable_as_double(_Tp __x) noexcept
822 {
823 static_assert(numeric_limits<_Tp>::is_integer, "");
824 static_assert(!numeric_limits<_Tp>::is_signed, "");
825 // All integers <= 2^53 are representable.
826 return (__x <= (1ull << __DBL_MANT_DIG__))
827 // Between 2^53 and 2^54 only even numbers are representable.
828 || (!(__x & 1) && __detail::__representable_as_double(__x >> 1));
829 }
830
831 // Determine whether x+1 is representable as double.
832 template<typename _Tp>
833 constexpr bool
834 __p1_representable_as_double(_Tp __x) noexcept
835 {
836 static_assert(numeric_limits<_Tp>::is_integer, "");
837 static_assert(!numeric_limits<_Tp>::is_signed, "");
838 return numeric_limits<_Tp>::digits < __DBL_MANT_DIG__
839 || (bool(__x + 1u) // return false if x+1 wraps around to zero
840 && __detail::__representable_as_double(__x + 1u));
841 }
842 }
843
844 template<typename _RandomNumberEngine, size_t __k>
848 {
849 constexpr result_type __range = max() - min();
850 size_t __j = __k;
851 const result_type __y = _M_y - min();
852#pragma GCC diagnostic push
853#pragma GCC diagnostic ignored "-Wc++17-extensions" // if constexpr
854 // Avoid using slower long double arithmetic if possible.
855 if constexpr (__detail::__p1_representable_as_double(__range))
856 __j *= __y / (__range + 1.0);
857 else
858 __j *= __y / (__range + 1.0L);
859#pragma GCC diagnostic pop
860 _M_y = _M_v[__j];
861 _M_v[__j] = _M_b();
862
863 return _M_y;
864 }
865
866 template<typename _RandomNumberEngine, size_t __k,
867 typename _CharT, typename _Traits>
871 {
872 using __ios_base = typename basic_ostream<_CharT, _Traits>::ios_base;
873
874 const typename __ios_base::fmtflags __flags = __os.flags();
875 const _CharT __fill = __os.fill();
876 const _CharT __space = __os.widen(' ');
877 __os.flags(__ios_base::dec | __ios_base::fixed | __ios_base::left);
878 __os.fill(__space);
879
880 __os << __x.base();
881 for (size_t __i = 0; __i < __k; ++__i)
882 __os << __space << __x._M_v[__i];
883 __os << __space << __x._M_y;
884
885 __os.flags(__flags);
886 __os.fill(__fill);
887 return __os;
888 }
889
890 template<typename _RandomNumberEngine, size_t __k,
891 typename _CharT, typename _Traits>
895 {
896 using __ios_base = typename basic_istream<_CharT, _Traits>::ios_base;
897
898 const typename __ios_base::fmtflags __flags = __is.flags();
899 __is.flags(__ios_base::dec | __ios_base::skipws);
900
901 __is >> __x._M_b;
902 for (size_t __i = 0; __i < __k; ++__i)
903 __is >> __x._M_v[__i];
904 __is >> __x._M_y;
905
906 __is.flags(__flags);
907 return __is;
908 }
909
910#if __glibcxx_philox_engine // >= C++26
911
912 template<typename _UIntType, size_t __w, size_t __n, size_t __r,
913 _UIntType... __consts>
914 _UIntType
915 philox_engine<_UIntType, __w, __n, __r, __consts...>::
916 _S_mulhi(_UIntType __a, _UIntType __b)
917 {
918 const __uint128_t __num =
919 static_cast<__uint128_t>(__a) * static_cast<__uint128_t>(__b);
920 return static_cast<_UIntType>((__num >> __w) & max());
921 }
922
923 template<typename _UIntType, size_t __w, size_t __n, size_t __r,
924 _UIntType... __consts>
925 _UIntType
926 philox_engine<_UIntType, __w, __n, __r, __consts...>::
927 _S_mullo(_UIntType __a, _UIntType __b)
928 {
929 return static_cast<_UIntType>((__a * __b) & max());
930 }
931
932 template<typename _UIntType, size_t __w, size_t __n, size_t __r,
933 _UIntType... __consts>
934 void
935 philox_engine<_UIntType, __w, __n, __r, __consts...>::_M_transition()
936 {
937 ++_M_i;
938 if (_M_i != __n)
939 return;
940
941 _M_philox();
942 if constexpr (__n == 4)
943 {
944 __uint128_t __uh =
945 (static_cast<__uint128_t>(_M_x[1]) << __w)
946 | (static_cast<__uint128_t>(_M_x[0]) + 1);
947 __uint128_t __lh =
948 ((static_cast<__uint128_t>(_M_x[3]) << __w)
949 | (_M_x[2]));
950 __uint128_t __bigMask =
951 (static_cast<__uint128_t>(1) << ((2 * __w) - 1))
952 | ((static_cast<__uint128_t>(1) << ((2 * __w) - 1)) - 1);
953 if ((__uh & __bigMask) == 0)
954 {
955 ++__lh;
956 __uh = 0;
957 }
958 _M_x[0] = __uh & max();
959 _M_x[1] = (__uh >> (__w)) & max();
960 _M_x[2] = __lh & max();
961 _M_x[3] = (__lh >> (__w)) & max();
962 }
963 else
964 {
965 __uint128_t __num =
966 (static_cast<__uint128_t>(_M_x[1]) << __w)
967 | (static_cast<__uint128_t>(_M_x[0]) + 1);
968 _M_x[0] = __num & max();
969 _M_x[1] = (__num >> __w) & max();
970 }
971 _M_i = 0;
972 }
973
974 template<typename _UIntType, size_t __w, size_t __n, size_t __r,
975 _UIntType... __consts>
976 void
977 philox_engine<_UIntType, __w, __n, __r, __consts...>::_M_philox()
978 {
979 array<_UIntType, __n> __outputSeq = _M_x;
980 for (size_t __j = 0; __j < __r; ++__j)
981 {
982 array<_UIntType, __n> __intermedSeq{};
983 if constexpr (__n == 4)
984 {
985 __intermedSeq[0] = __outputSeq[2];
986 __intermedSeq[1] = __outputSeq[1];
987 __intermedSeq[2] = __outputSeq[0];
988 __intermedSeq[3] = __outputSeq[3];
989 }
990 else
991 {
992 __intermedSeq[0] = __outputSeq[0];
993 __intermedSeq[1] = __outputSeq[1];
994 }
995 for (unsigned long __k = 0; __k < (__n/2); ++__k)
996 {
997 __outputSeq[2*__k]
998 = _S_mulhi(__intermedSeq[2*__k], multipliers[__k])
999 ^ (((_M_k[__k] + (__j * round_consts[__k])) & max()))
1000 ^ __intermedSeq[2*__k+1];
1001
1002 __outputSeq[(2*__k)+1]
1003 = _S_mullo(__intermedSeq[2*__k], multipliers[__k]);
1004 }
1005 }
1006 _M_y = __outputSeq;
1007 }
1008
1009 template<typename _UIntType, size_t __w, size_t __n, size_t __r,
1010 _UIntType... __consts>
1011 template<typename _Sseq>
1012 void
1013 philox_engine<_UIntType, __w, __n, __r, __consts...>::seed(_Sseq& __q)
1014 requires __is_seed_seq<_Sseq>
1015 {
1016 seed(0);
1017
1018 const unsigned __p = 1 + ((__w - 1) / 32);
1019 uint_least32_t __tmpArr[(__n/2) * __p];
1020 __q.generate(__tmpArr + 0, __tmpArr + ((__n/2) * __p));
1021 for (unsigned __k = 0; __k < (__n/2); ++__k)
1022 {
1023 unsigned long long __precalc = 0;
1024 for (unsigned __j = 0; __j < __p; ++__j)
1025 {
1026 unsigned long long __multiplicand = (1ull << (32 * __j));
1027 __precalc += (__tmpArr[__k * __p + __j] * __multiplicand) & max();
1028 }
1029 _M_k[__k] = __precalc;
1030 }
1031 }
1032#endif // philox_engine
1033
1034 template<typename _IntType, typename _CharT, typename _Traits>
1035 std::basic_ostream<_CharT, _Traits>&
1036 operator<<(std::basic_ostream<_CharT, _Traits>& __os,
1038 {
1039 using __ios_base = typename basic_ostream<_CharT, _Traits>::ios_base;
1040
1041 const typename __ios_base::fmtflags __flags = __os.flags();
1042 const _CharT __fill = __os.fill();
1043 const _CharT __space = __os.widen(' ');
1044 __os.flags(__ios_base::scientific | __ios_base::left);
1045 __os.fill(__space);
1046
1047 __os << __x.a() << __space << __x.b();
1048
1049 __os.flags(__flags);
1050 __os.fill(__fill);
1051 return __os;
1052 }
1053
1054 template<typename _IntType, typename _CharT, typename _Traits>
1055 std::basic_istream<_CharT, _Traits>&
1058 {
1059 using param_type
1061 using __ios_base = typename basic_istream<_CharT, _Traits>::ios_base;
1062
1063 const typename __ios_base::fmtflags __flags = __is.flags();
1064 __is.flags(__ios_base::dec | __ios_base::skipws);
1065
1066 _IntType __a, __b;
1067 if (__is >> __a >> __b)
1068 __x.param(param_type(__a, __b));
1069
1070 __is.flags(__flags);
1071 return __is;
1072 }
1073
1074
1075 template<typename _RealType>
1076 template<typename _ForwardIterator,
1077 typename _UniformRandomNumberGenerator>
1078 void
1080 __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
1081 _UniformRandomNumberGenerator& __urng,
1082 const param_type& __p)
1083 {
1084 __glibcxx_function_requires(_ForwardIteratorConcept<_ForwardIterator>)
1085 __detail::_Adaptor<_UniformRandomNumberGenerator, result_type>
1086 __aurng(__urng);
1087 auto __range = __p.b() - __p.a();
1088 while (__f != __t)
1089 *__f++ = __aurng() * __range + __p.a();
1090 }
1091
1092 template<typename _RealType, typename _CharT, typename _Traits>
1095 const uniform_real_distribution<_RealType>& __x)
1096 {
1097 using __ios_base = typename basic_ostream<_CharT, _Traits>::ios_base;
1098
1099 const typename __ios_base::fmtflags __flags = __os.flags();
1100 const _CharT __fill = __os.fill();
1101 const std::streamsize __precision = __os.precision();
1102 const _CharT __space = __os.widen(' ');
1103 __os.flags(__ios_base::scientific | __ios_base::left);
1104 __os.fill(__space);
1106
1107 __os << __x.a() << __space << __x.b();
1108
1109 __os.flags(__flags);
1110 __os.fill(__fill);
1111 __os.precision(__precision);
1112 return __os;
1113 }
1114
1115 template<typename _RealType, typename _CharT, typename _Traits>
1116 std::basic_istream<_CharT, _Traits>&
1119 {
1120 using param_type
1122 using __ios_base = typename basic_istream<_CharT, _Traits>::ios_base;
1123
1124 const typename __ios_base::fmtflags __flags = __is.flags();
1125 __is.flags(__ios_base::skipws);
1126
1127 _RealType __a, __b;
1128 if (__is >> __a >> __b)
1129 __x.param(param_type(__a, __b));
1130
1131 __is.flags(__flags);
1132 return __is;
1133 }
1134
1135
1136 template<typename _ForwardIterator,
1137 typename _UniformRandomNumberGenerator>
1138 void
1139 std::bernoulli_distribution::
1140 __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
1141 _UniformRandomNumberGenerator& __urng,
1142 const param_type& __p)
1143 {
1144 __glibcxx_function_requires(_ForwardIteratorConcept<_ForwardIterator>)
1145 __detail::_Adaptor<_UniformRandomNumberGenerator, double>
1146 __aurng(__urng);
1147 auto __limit = __p.p() * (__aurng.max() - __aurng.min());
1149 while (__f != __t)
1150 *__f++ = (__aurng() - __aurng.min()) < __limit;
1151 }
1152
1153 template<typename _CharT, typename _Traits>
1156 const bernoulli_distribution& __x)
1157 {
1158 using __ios_base = typename basic_ostream<_CharT, _Traits>::ios_base;
1159
1160 const typename __ios_base::fmtflags __flags = __os.flags();
1161 const _CharT __fill = __os.fill();
1162 const std::streamsize __precision = __os.precision();
1163 __os.flags(__ios_base::scientific | __ios_base::left);
1164 __os.fill(__os.widen(' '));
1166
1167 __os << __x.p();
1168
1169 __os.flags(__flags);
1170 __os.fill(__fill);
1171 __os.precision(__precision);
1172 return __os;
1173 }
1174
1175
1176 template<typename _IntType>
1177 template<typename _UniformRandomNumberGenerator>
1178 typename geometric_distribution<_IntType>::result_type
1180 operator()(_UniformRandomNumberGenerator& __urng,
1181 const param_type& __param)
1182 {
1183 // About the epsilon thing see this thread:
1184 // http://gcc.gnu.org/ml/gcc-patches/2006-10/msg00971.html
1185 const double __naf =
1187 // The largest _RealType convertible to _IntType.
1188 const double __thr =
1190 __detail::_Adaptor<_UniformRandomNumberGenerator, double>
1191 __aurng(__urng);
1192
1193 double __cand;
1194 do
1195 __cand = std::floor(std::log(1.0 - __aurng()) / __param._M_log_1_p);
1196 while (__cand >= __thr);
1197
1198 return result_type(__cand + __naf);
1199 }
1200
1201 template<typename _IntType>
1202 template<typename _ForwardIterator,
1203 typename _UniformRandomNumberGenerator>
1204 void
1206 __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
1207 _UniformRandomNumberGenerator& __urng,
1208 const param_type& __param)
1209 {
1210 __glibcxx_function_requires(_ForwardIteratorConcept<_ForwardIterator>)
1211 // About the epsilon thing see this thread:
1212 // http://gcc.gnu.org/ml/gcc-patches/2006-10/msg00971.html
1213 const double __naf =
1215 // The largest _RealType convertible to _IntType.
1216 const double __thr =
1218 __detail::_Adaptor<_UniformRandomNumberGenerator, double>
1219 __aurng(__urng);
1220
1221 while (__f != __t)
1222 {
1223 double __cand;
1224 do
1225 __cand = std::floor(std::log(1.0 - __aurng())
1226 / __param._M_log_1_p);
1227 while (__cand >= __thr);
1228
1229 *__f++ = __cand + __naf;
1230 }
1231 }
1232
1233 template<typename _IntType,
1234 typename _CharT, typename _Traits>
1237 const geometric_distribution<_IntType>& __x)
1238 {
1239 using __ios_base = typename basic_ostream<_CharT, _Traits>::ios_base;
1240
1241 const typename __ios_base::fmtflags __flags = __os.flags();
1242 const _CharT __fill = __os.fill();
1243 const std::streamsize __precision = __os.precision();
1244 __os.flags(__ios_base::scientific | __ios_base::left);
1245 __os.fill(__os.widen(' '));
1247
1248 __os << __x.p();
1249
1250 __os.flags(__flags);
1251 __os.fill(__fill);
1252 __os.precision(__precision);
1253 return __os;
1254 }
1255
1256 template<typename _IntType,
1257 typename _CharT, typename _Traits>
1258 std::basic_istream<_CharT, _Traits>&
1261 {
1262 using param_type = typename geometric_distribution<_IntType>::param_type;
1263 using __ios_base = typename basic_istream<_CharT, _Traits>::ios_base;
1264
1265 const typename __ios_base::fmtflags __flags = __is.flags();
1266 __is.flags(__ios_base::skipws);
1267
1268 double __p;
1269 if (__is >> __p)
1270 __x.param(param_type(__p));
1271
1272 __is.flags(__flags);
1273 return __is;
1274 }
1275
1276 // This is Leger's algorithm, also in Devroye, Ch. X, Example 1.5.
1277 template<typename _IntType>
1278 template<typename _UniformRandomNumberGenerator>
1279 typename negative_binomial_distribution<_IntType>::result_type
1281 operator()(_UniformRandomNumberGenerator& __urng)
1282 {
1283 const double __y = _M_gd(__urng);
1284
1285 // XXX Is the constructor too slow?
1287 return __poisson(__urng);
1288 }
1289
1290 template<typename _IntType>
1291 template<typename _UniformRandomNumberGenerator>
1294 operator()(_UniformRandomNumberGenerator& __urng,
1295 const param_type& __p)
1296 {
1298 param_type;
1299
1300 const double __y =
1301 _M_gd(__urng, param_type(__p.k(), (1.0 - __p.p()) / __p.p()));
1302
1304 return __poisson(__urng);
1305 }
1306
1307 template<typename _IntType>
1308 template<typename _ForwardIterator,
1309 typename _UniformRandomNumberGenerator>
1310 void
1311 negative_binomial_distribution<_IntType>::
1312 __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
1313 _UniformRandomNumberGenerator& __urng)
1314 {
1315 __glibcxx_function_requires(_ForwardIteratorConcept<_ForwardIterator>)
1316 while (__f != __t)
1317 {
1318 const double __y = _M_gd(__urng);
1319
1320 // XXX Is the constructor too slow?
1322 *__f++ = __poisson(__urng);
1323 }
1324 }
1325
1326 template<typename _IntType>
1327 template<typename _ForwardIterator,
1328 typename _UniformRandomNumberGenerator>
1329 void
1331 __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
1332 _UniformRandomNumberGenerator& __urng,
1333 const param_type& __p)
1334 {
1335 __glibcxx_function_requires(_ForwardIteratorConcept<_ForwardIterator>)
1336 typename std::gamma_distribution<result_type>::param_type
1337 __p2(__p.k(), (1.0 - __p.p()) / __p.p());
1338
1339 while (__f != __t)
1340 {
1341 const double __y = _M_gd(__urng, __p2);
1342
1343 std::poisson_distribution<result_type> __poisson(__y);
1344 *__f++ = __poisson(__urng);
1345 }
1346 }
1347
1348 template<typename _IntType, typename _CharT, typename _Traits>
1349 std::basic_ostream<_CharT, _Traits>&
1350 operator<<(std::basic_ostream<_CharT, _Traits>& __os,
1352 {
1353 using __ios_base = typename basic_ostream<_CharT, _Traits>::ios_base;
1354
1355 const typename __ios_base::fmtflags __flags = __os.flags();
1356 const _CharT __fill = __os.fill();
1357 const std::streamsize __precision = __os.precision();
1358 const _CharT __space = __os.widen(' ');
1359 __os.flags(__ios_base::scientific | __ios_base::left);
1360 __os.fill(__os.widen(' '));
1361 __os.precision(std::numeric_limits<double>::max_digits10);
1362
1363 __os << __x.k() << __space << __x.p()
1364 << __space << __x._M_gd;
1365
1366 __os.flags(__flags);
1367 __os.fill(__fill);
1368 __os.precision(__precision);
1369 return __os;
1370 }
1371
1372 template<typename _IntType, typename _CharT, typename _Traits>
1373 std::basic_istream<_CharT, _Traits>&
1374 operator>>(std::basic_istream<_CharT, _Traits>& __is,
1376 {
1377 using param_type
1379 using __ios_base = typename basic_istream<_CharT, _Traits>::ios_base;
1380
1381 const typename __ios_base::fmtflags __flags = __is.flags();
1382 __is.flags(__ios_base::skipws);
1383
1384 _IntType __k;
1385 double __p;
1386 if (__is >> __k >> __p >> __x._M_gd)
1387 __x.param(param_type(__k, __p));
1388
1389 __is.flags(__flags);
1390 return __is;
1391 }
1392
1393
1394 template<typename _IntType>
1395 void
1398 {
1399#if _GLIBCXX_USE_C99_MATH_FUNCS
1400 if (_M_mean >= 12)
1401 {
1402 const double __m = std::floor(_M_mean);
1403 _M_lm_thr = std::log(_M_mean);
1404 _M_lfm = std::lgamma(__m + 1);
1405 _M_sm = std::sqrt(__m);
1406
1407 const double __pi_4 = 0.7853981633974483096156608458198757L;
1408 const double __dx = std::sqrt(2 * __m * std::log(32 * __m
1409 / __pi_4));
1410 _M_d = std::round(std::max<double>(6.0, std::min(__m, __dx)));
1411 const double __cx = 2 * __m + _M_d;
1412 _M_scx = std::sqrt(__cx / 2);
1413 _M_1cx = 1 / __cx;
1414
1415 _M_c2b = std::sqrt(__pi_4 * __cx) * std::exp(_M_1cx);
1416 _M_cb = 2 * __cx * std::exp(-_M_d * _M_1cx * (1 + _M_d / 2))
1417 / _M_d;
1418 }
1419 else
1420#endif
1421 _M_lm_thr = std::exp(-_M_mean);
1422 }
1423
1424 /**
1425 * A rejection algorithm when mean >= 12 and a simple method based
1426 * upon the multiplication of uniform random variates otherwise.
1427 * NB: The former is available only if _GLIBCXX_USE_C99_MATH_FUNCS
1428 * is defined.
1429 *
1430 * Reference:
1431 * Devroye, L. Non-Uniform Random Variates Generation. Springer-Verlag,
1432 * New York, 1986, Ch. X, Sects. 3.3 & 3.4 (+ Errata!).
1433 */
1434 template<typename _IntType>
1435 template<typename _UniformRandomNumberGenerator>
1436 typename poisson_distribution<_IntType>::result_type
1438 operator()(_UniformRandomNumberGenerator& __urng,
1439 const param_type& __param)
1440 {
1441 __detail::_Adaptor<_UniformRandomNumberGenerator, double>
1442 __aurng(__urng);
1443#if _GLIBCXX_USE_C99_MATH_FUNCS
1444 if (__param.mean() >= 12)
1445 {
1446 double __x;
1447
1448 // See comments above...
1449 const double __naf =
1451 const double __thr =
1453
1454 const double __m = std::floor(__param.mean());
1455 // sqrt(pi / 2)
1456 const double __spi_2 = 1.2533141373155002512078826424055226L;
1457 const double __c1 = __param._M_sm * __spi_2;
1458 const double __c2 = __param._M_c2b + __c1;
1459 const double __c3 = __c2 + 1;
1460 const double __c4 = __c3 + 1;
1461 // 1 / 78
1462 const double __178 = 0.0128205128205128205128205128205128L;
1463 // e^(1 / 78)
1464 const double __e178 = 1.0129030479320018583185514777512983L;
1465 const double __c5 = __c4 + __e178;
1466 const double __c = __param._M_cb + __c5;
1467 const double __2cx = 2 * (2 * __m + __param._M_d);
1468
1469 bool __reject = true;
1470 do
1471 {
1472 const double __u = __c * __aurng();
1473 const double __e = -std::log(1.0 - __aurng());
1474
1475 double __w = 0.0;
1476
1477 if (__u <= __c1)
1478 {
1479 const double __n = _M_nd(__urng);
1480 const double __y = -std::abs(__n) * __param._M_sm - 1;
1481 __x = std::floor(__y);
1482 __w = -__n * __n / 2;
1483 if (__x < -__m)
1484 continue;
1485 }
1486 else if (__u <= __c2)
1487 {
1488 const double __n = _M_nd(__urng);
1489 const double __y = 1 + std::abs(__n) * __param._M_scx;
1490 __x = std::ceil(__y);
1491 __w = __y * (2 - __y) * __param._M_1cx;
1492 if (__x > __param._M_d)
1493 continue;
1494 }
1495 else if (__u <= __c3)
1496 // NB: This case not in the book, nor in the Errata,
1497 // but should be ok...
1498 __x = -1;
1499 else if (__u <= __c4)
1500 __x = 0;
1501 else if (__u <= __c5)
1502 {
1503 __x = 1;
1504 // Only in the Errata, see libstdc++/83237.
1505 __w = __178;
1506 }
1507 else
1508 {
1509 const double __v = -std::log(1.0 - __aurng());
1510 const double __y = __param._M_d
1511 + __v * __2cx / __param._M_d;
1512 __x = std::ceil(__y);
1513 __w = -__param._M_d * __param._M_1cx * (1 + __y / 2);
1514 }
1515
1516 __reject = (__w - __e - __x * __param._M_lm_thr
1517 > __param._M_lfm - std::lgamma(__x + __m + 1));
1518
1519 __reject |= __x + __m >= __thr;
1520
1521 } while (__reject);
1522
1523 return result_type(__x + __m + __naf);
1524 }
1525 else
1526#endif
1527 {
1528 _IntType __x = 0;
1529 double __prod = 1.0;
1530
1531 do
1532 {
1533 __prod *= __aurng();
1534 __x += 1;
1535 }
1536 while (__prod > __param._M_lm_thr);
1537
1538 return __x - 1;
1539 }
1540 }
1541
1542 template<typename _IntType>
1543 template<typename _ForwardIterator,
1544 typename _UniformRandomNumberGenerator>
1545 void
1547 __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
1548 _UniformRandomNumberGenerator& __urng,
1549 const param_type& __param)
1550 {
1551 __glibcxx_function_requires(_ForwardIteratorConcept<_ForwardIterator>)
1552 // We could duplicate everything from operator()...
1553 while (__f != __t)
1554 *__f++ = this->operator()(__urng, __param);
1555 }
1556
1557 template<typename _IntType,
1558 typename _CharT, typename _Traits>
1561 const poisson_distribution<_IntType>& __x)
1562 {
1563 using __ios_base = typename basic_ostream<_CharT, _Traits>::ios_base;
1564
1565 const typename __ios_base::fmtflags __flags = __os.flags();
1566 const _CharT __fill = __os.fill();
1567 const std::streamsize __precision = __os.precision();
1568 const _CharT __space = __os.widen(' ');
1569 __os.flags(__ios_base::scientific | __ios_base::left);
1570 __os.fill(__space);
1572
1573 __os << __x.mean() << __space << __x._M_nd;
1574
1575 __os.flags(__flags);
1576 __os.fill(__fill);
1577 __os.precision(__precision);
1578 return __os;
1579 }
1580
1581 template<typename _IntType,
1582 typename _CharT, typename _Traits>
1583 std::basic_istream<_CharT, _Traits>&
1584 operator>>(std::basic_istream<_CharT, _Traits>& __is,
1585 poisson_distribution<_IntType>& __x)
1586 {
1587 using param_type = typename poisson_distribution<_IntType>::param_type;
1588 using __ios_base = typename basic_istream<_CharT, _Traits>::ios_base;
1589
1590 const typename __ios_base::fmtflags __flags = __is.flags();
1591 __is.flags(__ios_base::skipws);
1592
1593 double __mean;
1594 if (__is >> __mean >> __x._M_nd)
1595 __x.param(param_type(__mean));
1596
1597 __is.flags(__flags);
1598 return __is;
1599 }
1600
1601
1602 template<typename _IntType>
1603 void
1606 {
1607 const double __p12 = _M_p <= 0.5 ? _M_p : 1.0 - _M_p;
1608
1609 _M_easy = true;
1610
1611#if _GLIBCXX_USE_C99_MATH_FUNCS
1612 if (_M_t * __p12 >= 8)
1613 {
1614 _M_easy = false;
1615 const double __np = std::floor(_M_t * __p12);
1616 const double __pa = __np / _M_t;
1617 const double __1p = 1 - __pa;
1618
1619 const double __pi_4 = 0.7853981633974483096156608458198757L;
1620 const double __d1x =
1621 std::sqrt(__np * __1p * std::log(32 * __np
1622 / (81 * __pi_4 * __1p)));
1623 _M_d1 = std::round(std::max<double>(1.0, __d1x));
1624 const double __d2x =
1625 std::sqrt(__np * __1p * std::log(32 * _M_t * __1p
1626 / (__pi_4 * __pa)));
1627 _M_d2 = std::round(std::max<double>(1.0, __d2x));
1628
1629 // sqrt(pi / 2)
1630 const double __spi_2 = 1.2533141373155002512078826424055226L;
1631 _M_s1 = std::sqrt(__np * __1p) * (1 + _M_d1 / (4 * __np));
1632 _M_s2 = std::sqrt(__np * __1p) * (1 + _M_d2 / (4 * (_M_t * __1p)));
1633 _M_c = 2 * _M_d1 / __np;
1634 _M_a1 = std::exp(_M_c) * _M_s1 * __spi_2;
1635 const double __a12 = _M_a1 + _M_s2 * __spi_2;
1636 const double __s1s = _M_s1 * _M_s1;
1637 _M_a123 = __a12 + (std::exp(_M_d1 / (_M_t * __1p))
1638 * 2 * __s1s / _M_d1
1639 * std::exp(-_M_d1 * _M_d1 / (2 * __s1s)));
1640 const double __s2s = _M_s2 * _M_s2;
1641 _M_s = (_M_a123 + 2 * __s2s / _M_d2
1642 * std::exp(-_M_d2 * _M_d2 / (2 * __s2s)));
1643 _M_lf = (std::lgamma(__np + 1)
1644 + std::lgamma(_M_t - __np + 1));
1645 _M_lp1p = std::log(__pa / __1p);
1646
1647 _M_q = -std::log(1 - (__p12 - __pa) / __1p);
1648 }
1649 else
1650#endif
1651 _M_q = -std::log(1 - __p12);
1652 }
1653
1654 template<typename _IntType>
1655 template<typename _UniformRandomNumberGenerator>
1658 _M_waiting(_UniformRandomNumberGenerator& __urng,
1659 _IntType __t, double __q)
1660 {
1661 _IntType __x = 0;
1662 double __sum = 0.0;
1663 __detail::_Adaptor<_UniformRandomNumberGenerator, double>
1664 __aurng(__urng);
1665
1666 do
1667 {
1668 if (__t == __x)
1669 return __x;
1670 const double __e = -std::log(1.0 - __aurng());
1671 __sum += __e / (__t - __x);
1672 __x += 1;
1673 }
1674 while (__sum <= __q);
1675
1676 return __x - 1;
1677 }
1678
1679 /**
1680 * A rejection algorithm when t * p >= 8 and a simple waiting time
1681 * method - the second in the referenced book - otherwise.
1682 * NB: The former is available only if _GLIBCXX_USE_C99_MATH_FUNCS
1683 * is defined.
1684 *
1685 * Reference:
1686 * Devroye, L. Non-Uniform Random Variates Generation. Springer-Verlag,
1687 * New York, 1986, Ch. X, Sect. 4 (+ Errata!).
1688 */
1689 template<typename _IntType>
1690 template<typename _UniformRandomNumberGenerator>
1693 operator()(_UniformRandomNumberGenerator& __urng,
1694 const param_type& __param)
1695 {
1696 result_type __ret;
1697 const _IntType __t = __param.t();
1698 const double __p = __param.p();
1699 const double __p12 = __p <= 0.5 ? __p : 1.0 - __p;
1700 __detail::_Adaptor<_UniformRandomNumberGenerator, double>
1701 __aurng(__urng);
1702
1703#if _GLIBCXX_USE_C99_MATH_FUNCS
1704 if (!__param._M_easy)
1705 {
1706 double __x;
1707
1708 // See comments above...
1709 const double __naf =
1711 const double __thr =
1713
1714 const double __np = std::floor(__t * __p12);
1715
1716 // sqrt(pi / 2)
1717 const double __spi_2 = 1.2533141373155002512078826424055226L;
1718 const double __a1 = __param._M_a1;
1719 const double __a12 = __a1 + __param._M_s2 * __spi_2;
1720 const double __a123 = __param._M_a123;
1721 const double __s1s = __param._M_s1 * __param._M_s1;
1722 const double __s2s = __param._M_s2 * __param._M_s2;
1723
1724 bool __reject;
1725 do
1726 {
1727 const double __u = __param._M_s * __aurng();
1728
1729 double __v;
1730
1731 if (__u <= __a1)
1732 {
1733 const double __n = _M_nd(__urng);
1734 const double __y = __param._M_s1 * std::abs(__n);
1735 __reject = __y >= __param._M_d1;
1736 if (!__reject)
1737 {
1738 const double __e = -std::log(1.0 - __aurng());
1739 __x = std::floor(__y);
1740 __v = -__e - __n * __n / 2 + __param._M_c;
1741 }
1742 }
1743 else if (__u <= __a12)
1744 {
1745 const double __n = _M_nd(__urng);
1746 const double __y = __param._M_s2 * std::abs(__n);
1747 __reject = __y >= __param._M_d2;
1748 if (!__reject)
1749 {
1750 const double __e = -std::log(1.0 - __aurng());
1751 __x = std::floor(-__y);
1752 __v = -__e - __n * __n / 2;
1753 }
1754 }
1755 else if (__u <= __a123)
1756 {
1757 const double __e1 = -std::log(1.0 - __aurng());
1758 const double __e2 = -std::log(1.0 - __aurng());
1759
1760 const double __y = __param._M_d1
1761 + 2 * __s1s * __e1 / __param._M_d1;
1762 __x = std::floor(__y);
1763 __v = (-__e2 + __param._M_d1 * (1 / (__t - __np)
1764 -__y / (2 * __s1s)));
1765 __reject = false;
1766 }
1767 else
1768 {
1769 const double __e1 = -std::log(1.0 - __aurng());
1770 const double __e2 = -std::log(1.0 - __aurng());
1771
1772 const double __y = __param._M_d2
1773 + 2 * __s2s * __e1 / __param._M_d2;
1774 __x = std::floor(-__y);
1775 __v = -__e2 - __param._M_d2 * __y / (2 * __s2s);
1776 __reject = false;
1777 }
1778
1779 __reject = __reject || __x < -__np || __x > __t - __np;
1780 if (!__reject)
1781 {
1782 const double __lfx =
1783 std::lgamma(__np + __x + 1)
1784 + std::lgamma(__t - (__np + __x) + 1);
1785 __reject = __v > __param._M_lf - __lfx
1786 + __x * __param._M_lp1p;
1787 }
1788
1789 __reject |= __x + __np >= __thr;
1790 }
1791 while (__reject);
1792
1793 __x += __np + __naf;
1794
1795 const _IntType __z = _M_waiting(__urng, __t - _IntType(__x),
1796 __param._M_q);
1797 __ret = _IntType(__x) + __z;
1798 }
1799 else
1800#endif
1801 __ret = _M_waiting(__urng, __t, __param._M_q);
1802
1803 if (__p12 != __p)
1804 __ret = __t - __ret;
1805 return __ret;
1806 }
1807
1808 template<typename _IntType>
1809 template<typename _ForwardIterator,
1810 typename _UniformRandomNumberGenerator>
1811 void
1813 __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
1814 _UniformRandomNumberGenerator& __urng,
1815 const param_type& __param)
1816 {
1817 __glibcxx_function_requires(_ForwardIteratorConcept<_ForwardIterator>)
1818 // We could duplicate everything from operator()...
1819 while (__f != __t)
1820 *__f++ = this->operator()(__urng, __param);
1821 }
1822
1823 template<typename _IntType,
1824 typename _CharT, typename _Traits>
1827 const binomial_distribution<_IntType>& __x)
1828 {
1829 using __ios_base = typename basic_ostream<_CharT, _Traits>::ios_base;
1830
1831 const typename __ios_base::fmtflags __flags = __os.flags();
1832 const _CharT __fill = __os.fill();
1833 const std::streamsize __precision = __os.precision();
1834 const _CharT __space = __os.widen(' ');
1835 __os.flags(__ios_base::scientific | __ios_base::left);
1836 __os.fill(__space);
1838
1839 __os << __x.t() << __space << __x.p()
1840 << __space << __x._M_nd;
1841
1842 __os.flags(__flags);
1843 __os.fill(__fill);
1844 __os.precision(__precision);
1845 return __os;
1846 }
1847
1848 template<typename _IntType,
1849 typename _CharT, typename _Traits>
1850 std::basic_istream<_CharT, _Traits>&
1851 operator>>(std::basic_istream<_CharT, _Traits>& __is,
1852 binomial_distribution<_IntType>& __x)
1853 {
1854 using param_type = typename binomial_distribution<_IntType>::param_type;
1855 using __ios_base = typename basic_istream<_CharT, _Traits>::ios_base;
1856
1857 const typename __ios_base::fmtflags __flags = __is.flags();
1858 __is.flags(__ios_base::dec | __ios_base::skipws);
1859
1860 _IntType __t;
1861 double __p;
1862 if (__is >> __t >> __p >> __x._M_nd)
1863 __x.param(param_type(__t, __p));
1864
1865 __is.flags(__flags);
1866 return __is;
1867 }
1868
1869
1870 template<typename _RealType>
1871 template<typename _ForwardIterator,
1872 typename _UniformRandomNumberGenerator>
1873 void
1875 __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
1876 _UniformRandomNumberGenerator& __urng,
1877 const param_type& __p)
1878 {
1879 __glibcxx_function_requires(_ForwardIteratorConcept<_ForwardIterator>)
1880 __detail::_Adaptor<_UniformRandomNumberGenerator, result_type>
1881 __aurng(__urng);
1882 while (__f != __t)
1883 *__f++ = -std::log(result_type(1) - __aurng()) / __p.lambda();
1884 }
1885
1886 template<typename _RealType, typename _CharT, typename _Traits>
1889 const exponential_distribution<_RealType>& __x)
1890 {
1891 using __ios_base = typename basic_ostream<_CharT, _Traits>::ios_base;
1892
1893 const typename __ios_base::fmtflags __flags = __os.flags();
1894 const _CharT __fill = __os.fill();
1895 const std::streamsize __precision = __os.precision();
1896 __os.flags(__ios_base::scientific | __ios_base::left);
1897 __os.fill(__os.widen(' '));
1899
1900 __os << __x.lambda();
1901
1902 __os.flags(__flags);
1903 __os.fill(__fill);
1904 __os.precision(__precision);
1905 return __os;
1906 }
1907
1908 template<typename _RealType, typename _CharT, typename _Traits>
1909 std::basic_istream<_CharT, _Traits>&
1912 {
1913 using param_type
1915 using __ios_base = typename basic_istream<_CharT, _Traits>::ios_base;
1916
1917 const typename __ios_base::fmtflags __flags = __is.flags();
1918 __is.flags(__ios_base::dec | __ios_base::skipws);
1919
1920 _RealType __lambda;
1921 if (__is >> __lambda)
1922 __x.param(param_type(__lambda));
1923
1924 __is.flags(__flags);
1925 return __is;
1926 }
1927
1928
1929 /**
1930 * Polar method due to Marsaglia.
1931 *
1932 * Devroye, L. Non-Uniform Random Variates Generation. Springer-Verlag,
1933 * New York, 1986, Ch. V, Sect. 4.4.
1934 */
1935 template<typename _RealType>
1936 template<typename _UniformRandomNumberGenerator>
1939 operator()(_UniformRandomNumberGenerator& __urng,
1940 const param_type& __param)
1941 {
1942 result_type __ret;
1943 __detail::_Adaptor<_UniformRandomNumberGenerator, result_type>
1944 __aurng(__urng);
1945
1946 if (_M_saved_available)
1947 {
1948 _M_saved_available = false;
1949 __ret = _M_saved;
1950 }
1951 else
1952 {
1953 result_type __x, __y, __r2;
1954 do
1955 {
1956 __x = result_type(2.0) * __aurng() - 1.0;
1957 __y = result_type(2.0) * __aurng() - 1.0;
1958 __r2 = __x * __x + __y * __y;
1959 }
1960 while (__r2 > 1.0 || __r2 == 0.0);
1961
1962 const result_type __mult = std::sqrt(-2 * std::log(__r2) / __r2);
1963 _M_saved = __x * __mult;
1964 _M_saved_available = true;
1965 __ret = __y * __mult;
1966 }
1967
1968 __ret = __ret * __param.stddev() + __param.mean();
1969 return __ret;
1970 }
1971
1972 template<typename _RealType>
1973 template<typename _ForwardIterator,
1974 typename _UniformRandomNumberGenerator>
1975 void
1977 __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
1978 _UniformRandomNumberGenerator& __urng,
1979 const param_type& __param)
1980 {
1981 __glibcxx_function_requires(_ForwardIteratorConcept<_ForwardIterator>)
1982
1983 if (__f == __t)
1984 return;
1985
1986 if (_M_saved_available)
1987 {
1988 _M_saved_available = false;
1989 *__f++ = _M_saved * __param.stddev() + __param.mean();
1990
1991 if (__f == __t)
1992 return;
1993 }
1994
1995 __detail::_Adaptor<_UniformRandomNumberGenerator, result_type>
1996 __aurng(__urng);
1997
1998 while (__f + 1 < __t)
1999 {
2000 result_type __x, __y, __r2;
2001 do
2002 {
2003 __x = result_type(2.0) * __aurng() - 1.0;
2004 __y = result_type(2.0) * __aurng() - 1.0;
2005 __r2 = __x * __x + __y * __y;
2006 }
2007 while (__r2 > 1.0 || __r2 == 0.0);
2008
2009 const result_type __mult = std::sqrt(-2 * std::log(__r2) / __r2);
2010 *__f++ = __y * __mult * __param.stddev() + __param.mean();
2011 *__f++ = __x * __mult * __param.stddev() + __param.mean();
2012 }
2013
2014 if (__f != __t)
2015 {
2016 result_type __x, __y, __r2;
2017 do
2018 {
2019 __x = result_type(2.0) * __aurng() - 1.0;
2020 __y = result_type(2.0) * __aurng() - 1.0;
2021 __r2 = __x * __x + __y * __y;
2022 }
2023 while (__r2 > 1.0 || __r2 == 0.0);
2024
2025 const result_type __mult = std::sqrt(-2 * std::log(__r2) / __r2);
2026 _M_saved = __x * __mult;
2027 _M_saved_available = true;
2028 *__f = __y * __mult * __param.stddev() + __param.mean();
2029 }
2030 }
2031
2032 template<typename _RealType>
2033 bool
2034 operator==(const std::normal_distribution<_RealType>& __d1,
2035 const std::normal_distribution<_RealType>& __d2)
2036 {
2037 if (__d1._M_param == __d2._M_param
2038 && __d1._M_saved_available == __d2._M_saved_available)
2039 return __d1._M_saved_available ? __d1._M_saved == __d2._M_saved : true;
2040 else
2041 return false;
2042 }
2043
2044 template<typename _RealType, typename _CharT, typename _Traits>
2045 std::basic_ostream<_CharT, _Traits>&
2046 operator<<(std::basic_ostream<_CharT, _Traits>& __os,
2047 const normal_distribution<_RealType>& __x)
2048 {
2049 using __ios_base = typename basic_ostream<_CharT, _Traits>::ios_base;
2050
2051 const typename __ios_base::fmtflags __flags = __os.flags();
2052 const _CharT __fill = __os.fill();
2053 const std::streamsize __precision = __os.precision();
2054 const _CharT __space = __os.widen(' ');
2055 __os.flags(__ios_base::scientific | __ios_base::left);
2056 __os.fill(__space);
2058
2059 __os << __x.mean() << __space << __x.stddev()
2060 << __space << __x._M_saved_available;
2061 if (__x._M_saved_available)
2062 __os << __space << __x._M_saved;
2063
2064 __os.flags(__flags);
2065 __os.fill(__fill);
2066 __os.precision(__precision);
2067 return __os;
2068 }
2069
2070 template<typename _RealType, typename _CharT, typename _Traits>
2071 std::basic_istream<_CharT, _Traits>&
2072 operator>>(std::basic_istream<_CharT, _Traits>& __is,
2073 normal_distribution<_RealType>& __x)
2074 {
2075 using param_type = typename normal_distribution<_RealType>::param_type;
2076 using __ios_base = typename basic_istream<_CharT, _Traits>::ios_base;
2077
2078 const typename __ios_base::fmtflags __flags = __is.flags();
2079 __is.flags(__ios_base::dec | __ios_base::skipws);
2080
2081 double __mean, __stddev;
2082 bool __saved_avail;
2083 if (__is >> __mean >> __stddev >> __saved_avail)
2084 {
2085 if (!__saved_avail || (__is >> __x._M_saved))
2086 {
2087 __x._M_saved_available = __saved_avail;
2088 __x.param(param_type(__mean, __stddev));
2089 }
2090 }
2091
2092 __is.flags(__flags);
2093 return __is;
2094 }
2095
2096
2097 template<typename _RealType>
2098 template<typename _ForwardIterator,
2099 typename _UniformRandomNumberGenerator>
2100 void
2101 lognormal_distribution<_RealType>::
2102 __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
2103 _UniformRandomNumberGenerator& __urng,
2104 const param_type& __p)
2105 {
2106 __glibcxx_function_requires(_ForwardIteratorConcept<_ForwardIterator>)
2107 while (__f != __t)
2108 *__f++ = std::exp(__p.s() * _M_nd(__urng) + __p.m());
2109 }
2110
2111 template<typename _RealType, typename _CharT, typename _Traits>
2112 std::basic_ostream<_CharT, _Traits>&
2113 operator<<(std::basic_ostream<_CharT, _Traits>& __os,
2114 const lognormal_distribution<_RealType>& __x)
2115 {
2116 using __ios_base = typename basic_ostream<_CharT, _Traits>::ios_base;
2117
2118 const typename __ios_base::fmtflags __flags = __os.flags();
2119 const _CharT __fill = __os.fill();
2120 const std::streamsize __precision = __os.precision();
2121 const _CharT __space = __os.widen(' ');
2122 __os.flags(__ios_base::scientific | __ios_base::left);
2123 __os.fill(__space);
2125
2126 __os << __x.m() << __space << __x.s()
2127 << __space << __x._M_nd;
2128
2129 __os.flags(__flags);
2130 __os.fill(__fill);
2131 __os.precision(__precision);
2132 return __os;
2133 }
2134
2135 template<typename _RealType, typename _CharT, typename _Traits>
2136 std::basic_istream<_CharT, _Traits>&
2137 operator>>(std::basic_istream<_CharT, _Traits>& __is,
2138 lognormal_distribution<_RealType>& __x)
2139 {
2140 using param_type
2141 = typename lognormal_distribution<_RealType>::param_type;
2142 using __ios_base = typename basic_istream<_CharT, _Traits>::ios_base;
2143
2144 const typename __ios_base::fmtflags __flags = __is.flags();
2145 __is.flags(__ios_base::dec | __ios_base::skipws);
2146
2147 _RealType __m, __s;
2148 if (__is >> __m >> __s >> __x._M_nd)
2149 __x.param(param_type(__m, __s));
2150
2151 __is.flags(__flags);
2152 return __is;
2153 }
2154
2155 template<typename _RealType>
2156 template<typename _ForwardIterator,
2157 typename _UniformRandomNumberGenerator>
2158 void
2159 std::chi_squared_distribution<_RealType>::
2160 __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
2161 _UniformRandomNumberGenerator& __urng)
2162 {
2163 __glibcxx_function_requires(_ForwardIteratorConcept<_ForwardIterator>)
2164 while (__f != __t)
2165 *__f++ = 2 * _M_gd(__urng);
2166 }
2167
2168 template<typename _RealType>
2169 template<typename _ForwardIterator,
2170 typename _UniformRandomNumberGenerator>
2171 void
2172 std::chi_squared_distribution<_RealType>::
2173 __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
2174 _UniformRandomNumberGenerator& __urng,
2175 const typename
2176 std::gamma_distribution<result_type>::param_type& __p)
2177 {
2178 __glibcxx_function_requires(_ForwardIteratorConcept<_ForwardIterator>)
2179 while (__f != __t)
2180 *__f++ = 2 * _M_gd(__urng, __p);
2181 }
2182
2183 template<typename _RealType, typename _CharT, typename _Traits>
2184 std::basic_ostream<_CharT, _Traits>&
2185 operator<<(std::basic_ostream<_CharT, _Traits>& __os,
2186 const chi_squared_distribution<_RealType>& __x)
2187 {
2188 using __ios_base = typename basic_ostream<_CharT, _Traits>::ios_base;
2189
2190 const typename __ios_base::fmtflags __flags = __os.flags();
2191 const _CharT __fill = __os.fill();
2192 const std::streamsize __precision = __os.precision();
2193 const _CharT __space = __os.widen(' ');
2194 __os.flags(__ios_base::scientific | __ios_base::left);
2195 __os.fill(__space);
2197
2198 __os << __x.n() << __space << __x._M_gd;
2199
2200 __os.flags(__flags);
2201 __os.fill(__fill);
2202 __os.precision(__precision);
2203 return __os;
2204 }
2205
2206 template<typename _RealType, typename _CharT, typename _Traits>
2207 std::basic_istream<_CharT, _Traits>&
2208 operator>>(std::basic_istream<_CharT, _Traits>& __is,
2209 chi_squared_distribution<_RealType>& __x)
2210 {
2211 using param_type
2212 = typename chi_squared_distribution<_RealType>::param_type;
2213 using __ios_base = typename basic_istream<_CharT, _Traits>::ios_base;
2214
2215 const typename __ios_base::fmtflags __flags = __is.flags();
2216 __is.flags(__ios_base::dec | __ios_base::skipws);
2217
2218 _RealType __n;
2219 if (__is >> __n >> __x._M_gd)
2220 __x.param(param_type(__n));
2221
2222 __is.flags(__flags);
2223 return __is;
2224 }
2225
2226
2227 template<typename _RealType>
2228 template<typename _UniformRandomNumberGenerator>
2231 operator()(_UniformRandomNumberGenerator& __urng,
2232 const param_type& __p)
2233 {
2234 __detail::_Adaptor<_UniformRandomNumberGenerator, result_type>
2235 __aurng(__urng);
2236 _RealType __u;
2237 do
2238 __u = __aurng();
2239 while (__u == 0.5);
2240
2241 const _RealType __pi = 3.1415926535897932384626433832795029L;
2242 return __p.a() + __p.b() * std::tan(__pi * __u);
2243 }
2244
2245 template<typename _RealType>
2246 template<typename _ForwardIterator,
2247 typename _UniformRandomNumberGenerator>
2248 void
2250 __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
2251 _UniformRandomNumberGenerator& __urng,
2252 const param_type& __p)
2253 {
2254 __glibcxx_function_requires(_ForwardIteratorConcept<_ForwardIterator>)
2255 const _RealType __pi = 3.1415926535897932384626433832795029L;
2256 __detail::_Adaptor<_UniformRandomNumberGenerator, result_type>
2257 __aurng(__urng);
2258 while (__f != __t)
2259 {
2260 _RealType __u;
2261 do
2262 __u = __aurng();
2263 while (__u == 0.5);
2264
2265 *__f++ = __p.a() + __p.b() * std::tan(__pi * __u);
2266 }
2267 }
2268
2269 template<typename _RealType, typename _CharT, typename _Traits>
2272 const cauchy_distribution<_RealType>& __x)
2273 {
2274 using __ios_base = typename basic_ostream<_CharT, _Traits>::ios_base;
2275
2276 const typename __ios_base::fmtflags __flags = __os.flags();
2277 const _CharT __fill = __os.fill();
2278 const std::streamsize __precision = __os.precision();
2279 const _CharT __space = __os.widen(' ');
2280 __os.flags(__ios_base::scientific | __ios_base::left);
2281 __os.fill(__space);
2283
2284 __os << __x.a() << __space << __x.b();
2285
2286 __os.flags(__flags);
2287 __os.fill(__fill);
2288 __os.precision(__precision);
2289 return __os;
2290 }
2291
2292 template<typename _RealType, typename _CharT, typename _Traits>
2293 std::basic_istream<_CharT, _Traits>&
2296 {
2297 using param_type = typename cauchy_distribution<_RealType>::param_type;
2298 using __ios_base = typename basic_istream<_CharT, _Traits>::ios_base;
2299
2300 const typename __ios_base::fmtflags __flags = __is.flags();
2301 __is.flags(__ios_base::dec | __ios_base::skipws);
2302
2303 _RealType __a, __b;
2304 if (__is >> __a >> __b)
2305 __x.param(param_type(__a, __b));
2306
2307 __is.flags(__flags);
2308 return __is;
2309 }
2310
2311
2312 template<typename _RealType>
2313 template<typename _ForwardIterator,
2314 typename _UniformRandomNumberGenerator>
2315 void
2317 __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
2318 _UniformRandomNumberGenerator& __urng)
2319 {
2320 __glibcxx_function_requires(_ForwardIteratorConcept<_ForwardIterator>)
2321 while (__f != __t)
2322 *__f++ = ((_M_gd_x(__urng) * n()) / (_M_gd_y(__urng) * m()));
2323 }
2324
2325 template<typename _RealType>
2326 template<typename _ForwardIterator,
2327 typename _UniformRandomNumberGenerator>
2328 void
2329 std::fisher_f_distribution<_RealType>::
2330 __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
2331 _UniformRandomNumberGenerator& __urng,
2332 const param_type& __p)
2333 {
2334 __glibcxx_function_requires(_ForwardIteratorConcept<_ForwardIterator>)
2335 typedef typename std::gamma_distribution<result_type>::param_type
2336 param_type;
2337 param_type __p1(__p.m() / 2);
2338 param_type __p2(__p.n() / 2);
2339 while (__f != __t)
2340 *__f++ = ((_M_gd_x(__urng, __p1) * n())
2341 / (_M_gd_y(__urng, __p2) * m()));
2342 }
2343
2344 template<typename _RealType, typename _CharT, typename _Traits>
2345 std::basic_ostream<_CharT, _Traits>&
2346 operator<<(std::basic_ostream<_CharT, _Traits>& __os,
2347 const fisher_f_distribution<_RealType>& __x)
2348 {
2349 using __ios_base = typename basic_ostream<_CharT, _Traits>::ios_base;
2350
2351 const typename __ios_base::fmtflags __flags = __os.flags();
2352 const _CharT __fill = __os.fill();
2353 const std::streamsize __precision = __os.precision();
2354 const _CharT __space = __os.widen(' ');
2355 __os.flags(__ios_base::scientific | __ios_base::left);
2356 __os.fill(__space);
2358
2359 __os << __x.m() << __space << __x.n()
2360 << __space << __x._M_gd_x << __space << __x._M_gd_y;
2361
2362 __os.flags(__flags);
2363 __os.fill(__fill);
2364 __os.precision(__precision);
2365 return __os;
2366 }
2367
2368 template<typename _RealType, typename _CharT, typename _Traits>
2369 std::basic_istream<_CharT, _Traits>&
2370 operator>>(std::basic_istream<_CharT, _Traits>& __is,
2371 fisher_f_distribution<_RealType>& __x)
2372 {
2373 using param_type
2374 = typename fisher_f_distribution<_RealType>::param_type;
2375 using __ios_base = typename basic_istream<_CharT, _Traits>::ios_base;
2376
2377 const typename __ios_base::fmtflags __flags = __is.flags();
2378 __is.flags(__ios_base::dec | __ios_base::skipws);
2379
2380 _RealType __m, __n;
2381 if (__is >> __m >> __n >> __x._M_gd_x >> __x._M_gd_y)
2382 __x.param(param_type(__m, __n));
2383
2384 __is.flags(__flags);
2385 return __is;
2386 }
2387
2388
2389 template<typename _RealType>
2390 template<typename _ForwardIterator,
2391 typename _UniformRandomNumberGenerator>
2392 void
2393 std::student_t_distribution<_RealType>::
2394 __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
2395 _UniformRandomNumberGenerator& __urng)
2396 {
2397 __glibcxx_function_requires(_ForwardIteratorConcept<_ForwardIterator>)
2398 while (__f != __t)
2399 *__f++ = _M_nd(__urng) * std::sqrt(n() / _M_gd(__urng));
2400 }
2401
2402 template<typename _RealType>
2403 template<typename _ForwardIterator,
2404 typename _UniformRandomNumberGenerator>
2405 void
2406 std::student_t_distribution<_RealType>::
2407 __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
2408 _UniformRandomNumberGenerator& __urng,
2409 const param_type& __p)
2410 {
2411 __glibcxx_function_requires(_ForwardIteratorConcept<_ForwardIterator>)
2412 typename std::gamma_distribution<result_type>::param_type
2413 __p2(__p.n() / 2, 2);
2414 while (__f != __t)
2415 *__f++ = _M_nd(__urng) * std::sqrt(__p.n() / _M_gd(__urng, __p2));
2416 }
2417
2418 template<typename _RealType, typename _CharT, typename _Traits>
2419 std::basic_ostream<_CharT, _Traits>&
2420 operator<<(std::basic_ostream<_CharT, _Traits>& __os,
2421 const student_t_distribution<_RealType>& __x)
2422 {
2423 using __ios_base = typename basic_ostream<_CharT, _Traits>::ios_base;
2424
2425 const typename __ios_base::fmtflags __flags = __os.flags();
2426 const _CharT __fill = __os.fill();
2427 const std::streamsize __precision = __os.precision();
2428 const _CharT __space = __os.widen(' ');
2429 __os.flags(__ios_base::scientific | __ios_base::left);
2430 __os.fill(__space);
2432
2433 __os << __x.n() << __space << __x._M_nd << __space << __x._M_gd;
2434
2435 __os.flags(__flags);
2436 __os.fill(__fill);
2437 __os.precision(__precision);
2438 return __os;
2439 }
2440
2441 template<typename _RealType, typename _CharT, typename _Traits>
2442 std::basic_istream<_CharT, _Traits>&
2443 operator>>(std::basic_istream<_CharT, _Traits>& __is,
2444 student_t_distribution<_RealType>& __x)
2445 {
2446 using param_type
2447 = typename student_t_distribution<_RealType>::param_type;
2448 using __ios_base = typename basic_istream<_CharT, _Traits>::ios_base;
2449
2450 const typename __ios_base::fmtflags __flags = __is.flags();
2451 __is.flags(__ios_base::dec | __ios_base::skipws);
2452
2453 _RealType __n;
2454 if (__is >> __n >> __x._M_nd >> __x._M_gd)
2455 __x.param(param_type(__n));
2456
2457 __is.flags(__flags);
2458 return __is;
2459 }
2460
2461
2462 template<typename _RealType>
2463 void
2464 gamma_distribution<_RealType>::param_type::
2465 _M_initialize()
2466 {
2467 _M_malpha = _M_alpha < 1.0 ? _M_alpha + _RealType(1.0) : _M_alpha;
2468
2469 const _RealType __a1 = _M_malpha - _RealType(1.0) / _RealType(3.0);
2470 _M_a2 = _RealType(1.0) / std::sqrt(_RealType(9.0) * __a1);
2471 }
2472
2473 /**
2474 * Marsaglia, G. and Tsang, W. W.
2475 * "A Simple Method for Generating Gamma Variables"
2476 * ACM Transactions on Mathematical Software, 26, 3, 363-372, 2000.
2477 */
2478 template<typename _RealType>
2479 template<typename _UniformRandomNumberGenerator>
2482 operator()(_UniformRandomNumberGenerator& __urng,
2483 const param_type& __param)
2484 {
2485 __detail::_Adaptor<_UniformRandomNumberGenerator, result_type>
2486 __aurng(__urng);
2487
2488 result_type __u, __v, __n;
2489 const result_type __a1 = (__param._M_malpha
2490 - _RealType(1.0) / _RealType(3.0));
2491
2492 do
2493 {
2494 do
2495 {
2496 __n = _M_nd(__urng);
2497 __v = result_type(1.0) + __param._M_a2 * __n;
2498 }
2499 while (__v <= 0.0);
2500
2501 __v = __v * __v * __v;
2502 __u = __aurng();
2503 }
2504 while (__u > result_type(1.0) - 0.0331 * __n * __n * __n * __n
2505 && (std::log(__u) > (0.5 * __n * __n + __a1
2506 * (1.0 - __v + std::log(__v)))));
2507
2508 if (__param.alpha() == __param._M_malpha)
2509 return __a1 * __v * __param.beta();
2510 else
2511 {
2512 do
2513 __u = __aurng();
2514 while (__u == 0.0);
2515
2516 return (std::pow(__u, result_type(1.0) / __param.alpha())
2517 * __a1 * __v * __param.beta());
2519 }
2520
2521 template<typename _RealType>
2522 template<typename _ForwardIterator,
2523 typename _UniformRandomNumberGenerator>
2524 void
2525 gamma_distribution<_RealType>::
2526 __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
2527 _UniformRandomNumberGenerator& __urng,
2528 const param_type& __param)
2529 {
2530 __glibcxx_function_requires(_ForwardIteratorConcept<_ForwardIterator>)
2531 __detail::_Adaptor<_UniformRandomNumberGenerator, result_type>
2532 __aurng(__urng);
2533
2534 result_type __u, __v, __n;
2535 const result_type __a1 = (__param._M_malpha
2536 - _RealType(1.0) / _RealType(3.0));
2537
2538 if (__param.alpha() == __param._M_malpha)
2539 while (__f != __t)
2540 {
2541 do
2542 {
2543 do
2544 {
2545 __n = _M_nd(__urng);
2546 __v = result_type(1.0) + __param._M_a2 * __n;
2547 }
2548 while (__v <= 0.0);
2549
2550 __v = __v * __v * __v;
2551 __u = __aurng();
2552 }
2553 while (__u > result_type(1.0) - 0.0331 * __n * __n * __n * __n
2554 && (std::log(__u) > (0.5 * __n * __n + __a1
2555 * (1.0 - __v + std::log(__v)))));
2556
2557 *__f++ = __a1 * __v * __param.beta();
2558 }
2559 else
2560 while (__f != __t)
2561 {
2562 do
2563 {
2564 do
2565 {
2566 __n = _M_nd(__urng);
2567 __v = result_type(1.0) + __param._M_a2 * __n;
2568 }
2569 while (__v <= 0.0);
2570
2571 __v = __v * __v * __v;
2572 __u = __aurng();
2573 }
2574 while (__u > result_type(1.0) - 0.0331 * __n * __n * __n * __n
2575 && (std::log(__u) > (0.5 * __n * __n + __a1
2576 * (1.0 - __v + std::log(__v)))));
2577
2578 do
2579 __u = __aurng();
2580 while (__u == 0.0);
2581
2582 *__f++ = (std::pow(__u, result_type(1.0) / __param.alpha())
2583 * __a1 * __v * __param.beta());
2584 }
2585 }
2586
2587 template<typename _RealType, typename _CharT, typename _Traits>
2588 std::basic_ostream<_CharT, _Traits>&
2589 operator<<(std::basic_ostream<_CharT, _Traits>& __os,
2590 const gamma_distribution<_RealType>& __x)
2591 {
2592 using __ios_base = typename basic_ostream<_CharT, _Traits>::ios_base;
2593
2594 const typename __ios_base::fmtflags __flags = __os.flags();
2595 const _CharT __fill = __os.fill();
2596 const std::streamsize __precision = __os.precision();
2597 const _CharT __space = __os.widen(' ');
2598 __os.flags(__ios_base::scientific | __ios_base::left);
2599 __os.fill(__space);
2601
2602 __os << __x.alpha() << __space << __x.beta()
2603 << __space << __x._M_nd;
2604
2605 __os.flags(__flags);
2606 __os.fill(__fill);
2607 __os.precision(__precision);
2608 return __os;
2609 }
2610
2611 template<typename _RealType, typename _CharT, typename _Traits>
2612 std::basic_istream<_CharT, _Traits>&
2613 operator>>(std::basic_istream<_CharT, _Traits>& __is,
2614 gamma_distribution<_RealType>& __x)
2615 {
2616 using param_type = typename gamma_distribution<_RealType>::param_type;
2617 using __ios_base = typename basic_istream<_CharT, _Traits>::ios_base;
2618
2619 const typename __ios_base::fmtflags __flags = __is.flags();
2620 __is.flags(__ios_base::dec | __ios_base::skipws);
2621
2622 _RealType __alpha_val, __beta_val;
2623 if (__is >> __alpha_val >> __beta_val >> __x._M_nd)
2624 __x.param(param_type(__alpha_val, __beta_val));
2625
2626 __is.flags(__flags);
2627 return __is;
2628 }
2629
2630
2631 template<typename _RealType>
2632 template<typename _UniformRandomNumberGenerator>
2635 operator()(_UniformRandomNumberGenerator& __urng,
2636 const param_type& __p)
2637 {
2638 __detail::_Adaptor<_UniformRandomNumberGenerator, result_type>
2639 __aurng(__urng);
2640 return __p.b() * std::pow(-std::log(result_type(1) - __aurng()),
2641 result_type(1) / __p.a());
2642 }
2643
2644 template<typename _RealType>
2645 template<typename _ForwardIterator,
2646 typename _UniformRandomNumberGenerator>
2647 void
2649 __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
2650 _UniformRandomNumberGenerator& __urng,
2651 const param_type& __p)
2652 {
2653 __glibcxx_function_requires(_ForwardIteratorConcept<_ForwardIterator>)
2654 __detail::_Adaptor<_UniformRandomNumberGenerator, result_type>
2655 __aurng(__urng);
2656 auto __inv_a = result_type(1) / __p.a();
2657
2658 while (__f != __t)
2659 *__f++ = __p.b() * std::pow(-std::log(result_type(1) - __aurng()),
2660 __inv_a);
2661 }
2662
2663 template<typename _RealType, typename _CharT, typename _Traits>
2666 const weibull_distribution<_RealType>& __x)
2667 {
2668 using __ios_base = typename basic_ostream<_CharT, _Traits>::ios_base;
2669
2670 const typename __ios_base::fmtflags __flags = __os.flags();
2671 const _CharT __fill = __os.fill();
2672 const std::streamsize __precision = __os.precision();
2673 const _CharT __space = __os.widen(' ');
2674 __os.flags(__ios_base::scientific | __ios_base::left);
2675 __os.fill(__space);
2677
2678 __os << __x.a() << __space << __x.b();
2679
2680 __os.flags(__flags);
2681 __os.fill(__fill);
2682 __os.precision(__precision);
2683 return __os;
2684 }
2685
2686 template<typename _RealType, typename _CharT, typename _Traits>
2687 std::basic_istream<_CharT, _Traits>&
2690 {
2691 using param_type = typename weibull_distribution<_RealType>::param_type;
2692 using __ios_base = typename basic_istream<_CharT, _Traits>::ios_base;
2693
2694 const typename __ios_base::fmtflags __flags = __is.flags();
2695 __is.flags(__ios_base::dec | __ios_base::skipws);
2696
2697 _RealType __a, __b;
2698 if (__is >> __a >> __b)
2699 __x.param(param_type(__a, __b));
2700
2701 __is.flags(__flags);
2702 return __is;
2703 }
2704
2705
2706 template<typename _RealType>
2707 template<typename _UniformRandomNumberGenerator>
2710 operator()(_UniformRandomNumberGenerator& __urng,
2711 const param_type& __p)
2712 {
2713 __detail::_Adaptor<_UniformRandomNumberGenerator, result_type>
2714 __aurng(__urng);
2715 return __p.a() - __p.b() * std::log(-std::log(result_type(1)
2716 - __aurng()));
2717 }
2718
2719 template<typename _RealType>
2720 template<typename _ForwardIterator,
2721 typename _UniformRandomNumberGenerator>
2722 void
2724 __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
2725 _UniformRandomNumberGenerator& __urng,
2726 const param_type& __p)
2727 {
2728 __glibcxx_function_requires(_ForwardIteratorConcept<_ForwardIterator>)
2729 __detail::_Adaptor<_UniformRandomNumberGenerator, result_type>
2730 __aurng(__urng);
2731
2732 while (__f != __t)
2733 *__f++ = __p.a() - __p.b() * std::log(-std::log(result_type(1)
2734 - __aurng()));
2735 }
2736
2737 template<typename _RealType, typename _CharT, typename _Traits>
2740 const extreme_value_distribution<_RealType>& __x)
2741 {
2742 using __ios_base = typename basic_ostream<_CharT, _Traits>::ios_base;
2743
2744 const typename __ios_base::fmtflags __flags = __os.flags();
2745 const _CharT __fill = __os.fill();
2746 const std::streamsize __precision = __os.precision();
2747 const _CharT __space = __os.widen(' ');
2748 __os.flags(__ios_base::scientific | __ios_base::left);
2749 __os.fill(__space);
2751
2752 __os << __x.a() << __space << __x.b();
2753
2754 __os.flags(__flags);
2755 __os.fill(__fill);
2756 __os.precision(__precision);
2757 return __os;
2758 }
2759
2760 template<typename _RealType, typename _CharT, typename _Traits>
2761 std::basic_istream<_CharT, _Traits>&
2764 {
2765 using param_type
2767 using __ios_base = typename basic_istream<_CharT, _Traits>::ios_base;
2768
2769 const typename __ios_base::fmtflags __flags = __is.flags();
2770 __is.flags(__ios_base::dec | __ios_base::skipws);
2771
2772 _RealType __a, __b;
2773 if (__is >> __a >> __b)
2774 __x.param(param_type(__a, __b));
2775
2776 __is.flags(__flags);
2777 return __is;
2778 }
2779
2780
2781 template<typename _IntType>
2782 void
2783 discrete_distribution<_IntType>::param_type::
2784 _M_initialize()
2785 {
2786 if (_M_prob.size() < 2)
2787 {
2788 _M_prob.clear();
2789 return;
2790 }
2791
2792 const double __sum = std::accumulate(_M_prob.begin(),
2793 _M_prob.end(), 0.0);
2794 __glibcxx_assert(__sum > 0);
2795 // Now normalize the probabilites.
2796 __detail::__normalize(_M_prob.begin(), _M_prob.end(), _M_prob.begin(),
2797 __sum);
2798 // Accumulate partial sums.
2799 _M_cp.reserve(_M_prob.size());
2800 std::partial_sum(_M_prob.begin(), _M_prob.end(),
2801 std::back_inserter(_M_cp));
2802 // Make sure the last cumulative probability is one.
2803 _M_cp[_M_cp.size() - 1] = 1.0;
2804 }
2805
2806 template<typename _IntType>
2807 template<typename _Func>
2808 discrete_distribution<_IntType>::param_type::
2809 param_type(size_t __nw, double __xmin, double __xmax, _Func __fw)
2810 : _M_prob(), _M_cp()
2811 {
2812 const size_t __n = __nw == 0 ? 1 : __nw;
2813 const double __delta = (__xmax - __xmin) / __n;
2814
2815 _M_prob.reserve(__n);
2816 for (size_t __k = 0; __k < __nw; ++__k)
2817 _M_prob.push_back(__fw(__xmin + __k * __delta + 0.5 * __delta));
2818
2819 _M_initialize();
2820 }
2821
2822 template<typename _IntType>
2823 template<typename _UniformRandomNumberGenerator>
2824 typename discrete_distribution<_IntType>::result_type
2825 discrete_distribution<_IntType>::
2826 operator()(_UniformRandomNumberGenerator& __urng,
2827 const param_type& __param)
2828 {
2829 if (__param._M_cp.empty())
2830 return result_type(0);
2831
2832 __detail::_Adaptor<_UniformRandomNumberGenerator, double>
2833 __aurng(__urng);
2834
2835 const double __p = __aurng();
2836 auto __pos = std::lower_bound(__param._M_cp.begin(),
2837 __param._M_cp.end(), __p);
2838
2839 return __pos - __param._M_cp.begin();
2840 }
2841
2842 template<typename _IntType>
2843 template<typename _ForwardIterator,
2844 typename _UniformRandomNumberGenerator>
2845 void
2846 discrete_distribution<_IntType>::
2847 __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
2848 _UniformRandomNumberGenerator& __urng,
2849 const param_type& __param)
2850 {
2851 __glibcxx_function_requires(_ForwardIteratorConcept<_ForwardIterator>)
2852
2853 if (__param._M_cp.empty())
2854 {
2855 while (__f != __t)
2856 *__f++ = result_type(0);
2857 return;
2858 }
2859
2860 __detail::_Adaptor<_UniformRandomNumberGenerator, double>
2861 __aurng(__urng);
2862
2863 while (__f != __t)
2864 {
2865 const double __p = __aurng();
2866 auto __pos = std::lower_bound(__param._M_cp.begin(),
2867 __param._M_cp.end(), __p);
2868
2869 *__f++ = __pos - __param._M_cp.begin();
2870 }
2871 }
2872
2873 template<typename _IntType, typename _CharT, typename _Traits>
2876 const discrete_distribution<_IntType>& __x)
2877 {
2878 using __ios_base = typename basic_ostream<_CharT, _Traits>::ios_base;
2879
2880 const typename __ios_base::fmtflags __flags = __os.flags();
2881 const _CharT __fill = __os.fill();
2882 const std::streamsize __precision = __os.precision();
2883 const _CharT __space = __os.widen(' ');
2884 __os.flags(__ios_base::scientific | __ios_base::left);
2885 __os.fill(__space);
2887
2888 std::vector<double> __prob = __x.probabilities();
2889 __os << __prob.size();
2890 for (auto __dit = __prob.begin(); __dit != __prob.end(); ++__dit)
2891 __os << __space << *__dit;
2892
2893 __os.flags(__flags);
2894 __os.fill(__fill);
2895 __os.precision(__precision);
2896 return __os;
2897 }
2898
2899namespace __detail
2900{
2901 template<typename _ValT, typename _CharT, typename _Traits>
2902 basic_istream<_CharT, _Traits>&
2903 __extract_params(basic_istream<_CharT, _Traits>& __is,
2904 vector<_ValT>& __vals, size_t __n)
2905 {
2906 __vals.reserve(__n);
2907 while (__n--)
2908 {
2909 _ValT __val;
2910 if (__is >> __val)
2911 __vals.push_back(__val);
2912 else
2913 break;
2914 }
2915 return __is;
2916 }
2917} // namespace __detail
2918
2919 template<typename _IntType, typename _CharT, typename _Traits>
2922 discrete_distribution<_IntType>& __x)
2923 {
2924 using __ios_base = typename basic_istream<_CharT, _Traits>::ios_base;
2925
2926 const typename __ios_base::fmtflags __flags = __is.flags();
2927 __is.flags(__ios_base::dec | __ios_base::skipws);
2928
2929 size_t __n;
2930 if (__is >> __n)
2931 {
2932 std::vector<double> __prob_vec;
2933 if (__detail::__extract_params(__is, __prob_vec, __n))
2934 __x.param({__prob_vec.begin(), __prob_vec.end()});
2935 }
2936
2937 __is.flags(__flags);
2938 return __is;
2939 }
2940
2941
2942 template<typename _RealType>
2943 void
2944 piecewise_constant_distribution<_RealType>::param_type::
2945 _M_initialize()
2946 {
2947 if (_M_int.size() < 2
2948 || (_M_int.size() == 2
2949 && _M_int[0] == _RealType(0)
2950 && _M_int[1] == _RealType(1)))
2951 {
2952 _M_int.clear();
2953 _M_den.clear();
2954 return;
2955 }
2956
2957 const double __sum = std::accumulate(_M_den.begin(),
2958 _M_den.end(), 0.0);
2959 __glibcxx_assert(__sum > 0);
2960
2961 __detail::__normalize(_M_den.begin(), _M_den.end(), _M_den.begin(),
2962 __sum);
2963
2964 _M_cp.reserve(_M_den.size());
2965 std::partial_sum(_M_den.begin(), _M_den.end(),
2966 std::back_inserter(_M_cp));
2967
2968 // Make sure the last cumulative probability is one.
2969 _M_cp[_M_cp.size() - 1] = 1.0;
2970
2971 for (size_t __k = 0; __k < _M_den.size(); ++__k)
2972 _M_den[__k] /= _M_int[__k + 1] - _M_int[__k];
2973 }
2974
2975 template<typename _RealType>
2976 template<typename _InputIteratorB, typename _InputIteratorW>
2977 piecewise_constant_distribution<_RealType>::param_type::
2978 param_type(_InputIteratorB __bbegin,
2979 _InputIteratorB __bend,
2980 _InputIteratorW __wbegin)
2981 : _M_int(), _M_den(), _M_cp()
2982 {
2983 if (__bbegin != __bend)
2984 {
2985 for (;;)
2986 {
2987 _M_int.push_back(*__bbegin);
2988 ++__bbegin;
2989 if (__bbegin == __bend)
2990 break;
2991
2992 _M_den.push_back(*__wbegin);
2993 ++__wbegin;
2994 }
2995 }
2996
2997 _M_initialize();
2998 }
2999
3000 template<typename _RealType>
3001 template<typename _Func>
3004 : _M_int(), _M_den(), _M_cp()
3005 {
3006 _M_int.reserve(__bl.size());
3007 for (auto __biter = __bl.begin(); __biter != __bl.end(); ++__biter)
3008 _M_int.push_back(*__biter);
3009
3010 _M_den.reserve(_M_int.size() - 1);
3011 for (size_t __k = 0; __k < _M_int.size() - 1; ++__k)
3012 _M_den.push_back(__fw(0.5 * (_M_int[__k + 1] + _M_int[__k])));
3013
3014 _M_initialize();
3015 }
3016
3017 template<typename _RealType>
3018 template<typename _Func>
3019 piecewise_constant_distribution<_RealType>::param_type::
3020 param_type(size_t __nw, _RealType __xmin, _RealType __xmax, _Func __fw)
3021 : _M_int(), _M_den(), _M_cp()
3022 {
3023 const size_t __n = __nw == 0 ? 1 : __nw;
3024 const _RealType __delta = (__xmax - __xmin) / __n;
3025
3026 _M_int.reserve(__n + 1);
3027 for (size_t __k = 0; __k <= __nw; ++__k)
3028 _M_int.push_back(__xmin + __k * __delta);
3029
3030 _M_den.reserve(__n);
3031 for (size_t __k = 0; __k < __nw; ++__k)
3032 _M_den.push_back(__fw(_M_int[__k] + 0.5 * __delta));
3033
3034 _M_initialize();
3035 }
3036
3037 template<typename _RealType>
3038 template<typename _UniformRandomNumberGenerator>
3041 operator()(_UniformRandomNumberGenerator& __urng,
3042 const param_type& __param)
3043 {
3044 __detail::_Adaptor<_UniformRandomNumberGenerator, double>
3045 __aurng(__urng);
3046
3047 const double __p = __aurng();
3048 if (__param._M_cp.empty())
3049 return __p;
3050
3051 auto __pos = std::lower_bound(__param._M_cp.begin(),
3052 __param._M_cp.end(), __p);
3053 const size_t __i = __pos - __param._M_cp.begin();
3054
3055 const double __pref = __i > 0 ? __param._M_cp[__i - 1] : 0.0;
3056
3057 return __param._M_int[__i] + (__p - __pref) / __param._M_den[__i];
3058 }
3059
3060 template<typename _RealType>
3061 template<typename _ForwardIterator,
3062 typename _UniformRandomNumberGenerator>
3063 void
3065 __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
3066 _UniformRandomNumberGenerator& __urng,
3067 const param_type& __param)
3068 {
3069 __glibcxx_function_requires(_ForwardIteratorConcept<_ForwardIterator>)
3070 __detail::_Adaptor<_UniformRandomNumberGenerator, double>
3071 __aurng(__urng);
3072
3073 if (__param._M_cp.empty())
3074 {
3075 while (__f != __t)
3076 *__f++ = __aurng();
3077 return;
3078 }
3079
3080 while (__f != __t)
3081 {
3082 const double __p = __aurng();
3083
3084 auto __pos = std::lower_bound(__param._M_cp.begin(),
3085 __param._M_cp.end(), __p);
3086 const size_t __i = __pos - __param._M_cp.begin();
3087
3088 const double __pref = __i > 0 ? __param._M_cp[__i - 1] : 0.0;
3089
3090 *__f++ = (__param._M_int[__i]
3091 + (__p - __pref) / __param._M_den[__i]);
3092 }
3093 }
3094
3095 template<typename _RealType, typename _CharT, typename _Traits>
3096 std::basic_ostream<_CharT, _Traits>&
3097 operator<<(std::basic_ostream<_CharT, _Traits>& __os,
3099 {
3100 using __ios_base = typename basic_ostream<_CharT, _Traits>::ios_base;
3101
3102 const typename __ios_base::fmtflags __flags = __os.flags();
3103 const _CharT __fill = __os.fill();
3104 const std::streamsize __precision = __os.precision();
3105 const _CharT __space = __os.widen(' ');
3106 __os.flags(__ios_base::scientific | __ios_base::left);
3107 __os.fill(__space);
3109
3110 std::vector<_RealType> __int = __x.intervals();
3111 __os << __int.size() - 1;
3112
3113 for (auto __xit = __int.begin(); __xit != __int.end(); ++__xit)
3114 __os << __space << *__xit;
3115
3116 std::vector<double> __den = __x.densities();
3117 for (auto __dit = __den.begin(); __dit != __den.end(); ++__dit)
3118 __os << __space << *__dit;
3119
3120 __os.flags(__flags);
3121 __os.fill(__fill);
3122 __os.precision(__precision);
3123 return __os;
3124 }
3125
3126 template<typename _RealType, typename _CharT, typename _Traits>
3127 std::basic_istream<_CharT, _Traits>&
3128 operator>>(std::basic_istream<_CharT, _Traits>& __is,
3130 {
3131 using __ios_base = typename basic_istream<_CharT, _Traits>::ios_base;
3132
3133 const typename __ios_base::fmtflags __flags = __is.flags();
3134 __is.flags(__ios_base::dec | __ios_base::skipws);
3135
3136 size_t __n;
3137 if (__is >> __n)
3138 {
3139 std::vector<_RealType> __int_vec;
3140 if (__detail::__extract_params(__is, __int_vec, __n + 1))
3141 {
3142 std::vector<double> __den_vec;
3143 if (__detail::__extract_params(__is, __den_vec, __n))
3144 {
3145 __x.param({ __int_vec.begin(), __int_vec.end(),
3146 __den_vec.begin() });
3147 }
3148 }
3149 }
3150
3151 __is.flags(__flags);
3152 return __is;
3153 }
3154
3155
3156 template<typename _RealType>
3157 void
3160 {
3161 if (_M_int.size() < 2
3162 || (_M_int.size() == 2
3163 && _M_int[0] == _RealType(0)
3164 && _M_int[1] == _RealType(1)
3165 && _M_den[0] == _M_den[1]))
3166 {
3167 _M_int.clear();
3168 _M_den.clear();
3169 return;
3170 }
3171
3172 double __sum = 0.0;
3173 _M_cp.reserve(_M_int.size() - 1);
3174 _M_m.reserve(_M_int.size() - 1);
3175 for (size_t __k = 0; __k < _M_int.size() - 1; ++__k)
3176 {
3177 const _RealType __delta = _M_int[__k + 1] - _M_int[__k];
3178 __sum += 0.5 * (_M_den[__k + 1] + _M_den[__k]) * __delta;
3179 _M_cp.push_back(__sum);
3180 _M_m.push_back((_M_den[__k + 1] - _M_den[__k]) / __delta);
3181 }
3182 __glibcxx_assert(__sum > 0);
3183
3184 // Now normalize the densities...
3185 __detail::__normalize(_M_den.begin(), _M_den.end(), _M_den.begin(),
3186 __sum);
3187 // ... and partial sums...
3188 __detail::__normalize(_M_cp.begin(), _M_cp.end(), _M_cp.begin(), __sum);
3189 // ... and slopes.
3190 __detail::__normalize(_M_m.begin(), _M_m.end(), _M_m.begin(), __sum);
3191
3192 // Make sure the last cumulative probablility is one.
3193 _M_cp[_M_cp.size() - 1] = 1.0;
3194 }
3195
3196 template<typename _RealType>
3197 template<typename _InputIteratorB, typename _InputIteratorW>
3198 piecewise_linear_distribution<_RealType>::param_type::
3199 param_type(_InputIteratorB __bbegin,
3200 _InputIteratorB __bend,
3201 _InputIteratorW __wbegin)
3202 : _M_int(), _M_den(), _M_cp(), _M_m()
3203 {
3204 for (; __bbegin != __bend; ++__bbegin, (void) ++__wbegin)
3205 {
3206 _M_int.push_back(*__bbegin);
3207 _M_den.push_back(*__wbegin);
3208 }
3209
3210 _M_initialize();
3211 }
3212
3213 template<typename _RealType>
3214 template<typename _Func>
3217 : _M_int(), _M_den(), _M_cp(), _M_m()
3218 {
3219 _M_int.reserve(__bl.size());
3220 _M_den.reserve(__bl.size());
3221 for (auto __biter = __bl.begin(); __biter != __bl.end(); ++__biter)
3222 {
3223 _M_int.push_back(*__biter);
3224 _M_den.push_back(__fw(*__biter));
3225 }
3226
3227 _M_initialize();
3228 }
3229
3230 template<typename _RealType>
3231 template<typename _Func>
3233 param_type(size_t __nw, _RealType __xmin, _RealType __xmax, _Func __fw)
3234 : _M_int(), _M_den(), _M_cp(), _M_m()
3235 {
3236 const size_t __n = __nw == 0 ? 1 : __nw;
3237 const _RealType __delta = (__xmax - __xmin) / __n;
3238
3239 _M_int.reserve(__n + 1);
3240 _M_den.reserve(__n + 1);
3241 for (size_t __k = 0; __k <= __nw; ++__k)
3242 {
3243 _M_int.push_back(__xmin + __k * __delta);
3244 _M_den.push_back(__fw(_M_int[__k] + __delta));
3245 }
3246
3247 _M_initialize();
3248 }
3249
3250 template<typename _RealType>
3251 template<typename _UniformRandomNumberGenerator>
3254 operator()(_UniformRandomNumberGenerator& __urng,
3255 const param_type& __param)
3256 {
3257 __detail::_Adaptor<_UniformRandomNumberGenerator, double>
3258 __aurng(__urng);
3259
3260 const double __p = __aurng();
3261 if (__param._M_cp.empty())
3262 return __p;
3263
3264 auto __pos = std::lower_bound(__param._M_cp.begin(),
3265 __param._M_cp.end(), __p);
3266 const size_t __i = __pos - __param._M_cp.begin();
3267
3268 const double __pref = __i > 0 ? __param._M_cp[__i - 1] : 0.0;
3269
3270 const double __a = 0.5 * __param._M_m[__i];
3271 const double __b = __param._M_den[__i];
3272 const double __cm = __p - __pref;
3273
3274 _RealType __x = __param._M_int[__i];
3275 if (__a == 0)
3276 __x += __cm / __b;
3277 else
3278 {
3279 const double __d = __b * __b + 4.0 * __a * __cm;
3280 __x += 0.5 * (std::sqrt(__d) - __b) / __a;
3281 }
3282
3283 return __x;
3284 }
3285
3286 template<typename _RealType>
3287 template<typename _ForwardIterator,
3288 typename _UniformRandomNumberGenerator>
3289 void
3291 __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
3292 _UniformRandomNumberGenerator& __urng,
3293 const param_type& __param)
3294 {
3295 __glibcxx_function_requires(_ForwardIteratorConcept<_ForwardIterator>)
3296 // We could duplicate everything from operator()...
3297 while (__f != __t)
3298 *__f++ = this->operator()(__urng, __param);
3299 }
3300
3301 template<typename _RealType, typename _CharT, typename _Traits>
3302 std::basic_ostream<_CharT, _Traits>&
3303 operator<<(std::basic_ostream<_CharT, _Traits>& __os,
3305 {
3306 using __ios_base = typename basic_ostream<_CharT, _Traits>::ios_base;
3307
3308 const typename __ios_base::fmtflags __flags = __os.flags();
3309 const _CharT __fill = __os.fill();
3310 const std::streamsize __precision = __os.precision();
3311 const _CharT __space = __os.widen(' ');
3312 __os.flags(__ios_base::scientific | __ios_base::left);
3313 __os.fill(__space);
3315
3316 std::vector<_RealType> __int = __x.intervals();
3317 __os << __int.size() - 1;
3318
3319 for (auto __xit = __int.begin(); __xit != __int.end(); ++__xit)
3320 __os << __space << *__xit;
3321
3322 std::vector<double> __den = __x.densities();
3323 for (auto __dit = __den.begin(); __dit != __den.end(); ++__dit)
3324 __os << __space << *__dit;
3325
3326 __os.flags(__flags);
3327 __os.fill(__fill);
3328 __os.precision(__precision);
3329 return __os;
3330 }
3331
3332 template<typename _RealType, typename _CharT, typename _Traits>
3333 std::basic_istream<_CharT, _Traits>&
3334 operator>>(std::basic_istream<_CharT, _Traits>& __is,
3336 {
3337 using __ios_base = typename basic_istream<_CharT, _Traits>::ios_base;
3338
3339 const typename __ios_base::fmtflags __flags = __is.flags();
3340 __is.flags(__ios_base::dec | __ios_base::skipws);
3341
3342 size_t __n;
3343 if (__is >> __n)
3344 {
3345 vector<_RealType> __int_vec;
3346 if (__detail::__extract_params(__is, __int_vec, __n + 1))
3347 {
3348 vector<double> __den_vec;
3349 if (__detail::__extract_params(__is, __den_vec, __n + 1))
3350 {
3351 __x.param({ __int_vec.begin(), __int_vec.end(),
3352 __den_vec.begin() });
3353 }
3354 }
3355 }
3356 __is.flags(__flags);
3357 return __is;
3358 }
3359
3360
3361 template<typename _IntType, typename>
3362 seed_seq::seed_seq(std::initializer_list<_IntType> __il)
3363 {
3364 _M_v.reserve(__il.size());
3365 for (auto __iter = __il.begin(); __iter != __il.end(); ++__iter)
3366 _M_v.push_back(__detail::__mod<result_type,
3367 __detail::_Shift<result_type, 32>::__value>(*__iter));
3368 }
3369
3370 template<typename _InputIterator>
3371 seed_seq::seed_seq(_InputIterator __begin, _InputIterator __end)
3372 {
3373#pragma GCC diagnostic push
3374#pragma GCC diagnostic ignored "-Wc++17-extensions" // if constexpr
3375 if constexpr (__is_random_access_iter<_InputIterator>::value)
3376 _M_v.reserve(std::distance(__begin, __end));
3377#pragma GCC diagnostic pop
3378
3379 for (_InputIterator __iter = __begin; __iter != __end; ++__iter)
3380 _M_v.push_back(__detail::__mod<result_type,
3381 __detail::_Shift<result_type, 32>::__value>(*__iter));
3382 }
3383
3384 template<typename _RandomAccessIterator>
3385 void
3386 seed_seq::generate(_RandomAccessIterator __begin,
3387 _RandomAccessIterator __end)
3388 {
3389 typedef typename iterator_traits<_RandomAccessIterator>::value_type
3390 _Type;
3391
3392 if (__begin == __end)
3393 return;
3394
3395 std::fill(__begin, __end, _Type(0x8b8b8b8bu));
3396
3397 const size_t __n = __end - __begin;
3398 const size_t __s = _M_v.size();
3399 const size_t __t = (__n >= 623) ? 11
3400 : (__n >= 68) ? 7
3401 : (__n >= 39) ? 5
3402 : (__n >= 7) ? 3
3403 : (__n - 1) / 2;
3404 const size_t __p = (__n - __t) / 2;
3405 const size_t __q = __p + __t;
3406 const size_t __m = std::max(size_t(__s + 1), __n);
3407
3408#ifndef __UINT32_TYPE__
3409 struct _Up
3410 {
3411 _Up(uint_least32_t v) : _M_v(v & 0xffffffffu) { }
3412
3413 operator uint_least32_t() const { return _M_v; }
3414
3415 uint_least32_t _M_v;
3416 };
3417 using uint32_t = _Up;
3418#endif
3419
3420 // k == 0, every element in [begin,end) equals 0x8b8b8b8bu
3421 {
3422 uint32_t __r1 = 1371501266u;
3423 uint32_t __r2 = __r1 + __s;
3424 __begin[__p] += __r1;
3425 __begin[__q] = (uint32_t)__begin[__q] + __r2;
3426 __begin[0] = __r2;
3427 }
3428
3429 for (size_t __k = 1; __k <= __s; ++__k)
3430 {
3431 const size_t __kn = __k % __n;
3432 const size_t __kpn = (__k + __p) % __n;
3433 const size_t __kqn = (__k + __q) % __n;
3434 uint32_t __arg = (__begin[__kn]
3435 ^ __begin[__kpn]
3436 ^ __begin[(__k - 1) % __n]);
3437 uint32_t __r1 = 1664525u * (__arg ^ (__arg >> 27));
3438 uint32_t __r2 = __r1 + (uint32_t)__kn + _M_v[__k - 1];
3439 __begin[__kpn] = (uint32_t)__begin[__kpn] + __r1;
3440 __begin[__kqn] = (uint32_t)__begin[__kqn] + __r2;
3441 __begin[__kn] = __r2;
3442 }
3443
3444 for (size_t __k = __s + 1; __k < __m; ++__k)
3445 {
3446 const size_t __kn = __k % __n;
3447 const size_t __kpn = (__k + __p) % __n;
3448 const size_t __kqn = (__k + __q) % __n;
3449 uint32_t __arg = (__begin[__kn]
3450 ^ __begin[__kpn]
3451 ^ __begin[(__k - 1) % __n]);
3452 uint32_t __r1 = 1664525u * (__arg ^ (__arg >> 27));
3453 uint32_t __r2 = __r1 + (uint32_t)__kn;
3454 __begin[__kpn] = (uint32_t)__begin[__kpn] + __r1;
3455 __begin[__kqn] = (uint32_t)__begin[__kqn] + __r2;
3456 __begin[__kn] = __r2;
3457 }
3458
3459 for (size_t __k = __m; __k < __m + __n; ++__k)
3460 {
3461 const size_t __kn = __k % __n;
3462 const size_t __kpn = (__k + __p) % __n;
3463 const size_t __kqn = (__k + __q) % __n;
3464 uint32_t __arg = (__begin[__kn]
3465 + __begin[__kpn]
3466 + __begin[(__k - 1) % __n]);
3467 uint32_t __r3 = 1566083941u * (__arg ^ (__arg >> 27));
3468 uint32_t __r4 = __r3 - __kn;
3469 __begin[__kpn] ^= __r3;
3470 __begin[__kqn] ^= __r4;
3471 __begin[__kn] = __r4;
3472 }
3473 }
3474
3475 template<typename _RealType, size_t __bits,
3476 typename _UniformRandomNumberGenerator>
3477 _RealType
3478 generate_canonical(_UniformRandomNumberGenerator& __urng)
3479 {
3481 "template argument must be a floating point type");
3482
3483 const size_t __b
3484 = std::min(static_cast<size_t>(std::numeric_limits<_RealType>::digits),
3485 __bits);
3486 const long double __r = static_cast<long double>(__urng.max())
3487 - static_cast<long double>(__urng.min()) + 1.0L;
3488 const size_t __log2r = std::log(__r) / std::log(2.0L);
3489 const size_t __m = std::max<size_t>(1UL,
3490 (__b + __log2r - 1UL) / __log2r);
3491 _RealType __ret;
3492 _RealType __sum = _RealType(0);
3493 _RealType __tmp = _RealType(1);
3494 for (size_t __k = __m; __k != 0; --__k)
3495 {
3496 __sum += _RealType(__urng() - __urng.min()) * __tmp;
3497 __tmp *= __r;
3498 }
3499 __ret = __sum / __tmp;
3500 if (__builtin_expect(__ret >= _RealType(1), 0))
3501 {
3502#if _GLIBCXX_USE_C99_MATH_FUNCS
3503 __ret = std::nextafter(_RealType(1), _RealType(0));
3504#else
3505 __ret = _RealType(1)
3506 - std::numeric_limits<_RealType>::epsilon() / _RealType(2);
3507#endif
3508 }
3509 return __ret;
3510 }
3511
3512_GLIBCXX_END_NAMESPACE_VERSION
3513} // namespace
3514
3515#endif
complex< _Tp > log(const complex< _Tp > &)
Return complex natural logarithm of z.
Definition complex:1162
complex< _Tp > tan(const complex< _Tp > &)
Return complex tangent of z.
Definition complex:1298
_Tp abs(const complex< _Tp > &)
Return magnitude of z.
Definition complex:968
complex< _Tp > exp(const complex< _Tp > &)
Return complex base e exponential of z.
Definition complex:1135
complex< _Tp > pow(const complex< _Tp > &, int)
Return x to the y'th power.
Definition complex:1357
complex< _Tp > sqrt(const complex< _Tp > &)
Return complex square root of z.
Definition complex:1271
constexpr const _Tp & max(const _Tp &, const _Tp &)
This does what you think it does.
constexpr const _Tp & min(const _Tp &, const _Tp &)
This does what you think it does.
_RealType generate_canonical(_UniformRandomNumberGenerator &__g)
A function template for converting the output of a (integral) uniform random number generator to a fl...
constexpr back_insert_iterator< _Container > back_inserter(_Container &__x)
constexpr _Tp accumulate(_InputIterator __first, _InputIterator __last, _Tp __init)
Accumulate values in a range.
constexpr _OutputIterator partial_sum(_InputIterator __first, _InputIterator __last, _OutputIterator __result)
Return list of partial sums.
ISO C++ entities toplevel namespace is std.
ptrdiff_t streamsize
Integral type for I/O operation counts and buffer sizes.
Definition postypes.h:73
constexpr iterator_traits< _InputIterator >::difference_type distance(_InputIterator __first, _InputIterator __last)
A generalization of pointer arithmetic.
constexpr _Tp __lg(_Tp __n)
This is a helper function for the sort routines and for random.tcc.
std::basic_istream< _CharT, _Traits > & operator>>(std::basic_istream< _CharT, _Traits > &__is, bitset< _Nb > &__x)
Global I/O operators for bitsets.
Definition bitset:1637
std::basic_ostream< _CharT, _Traits > & operator<<(std::basic_ostream< _CharT, _Traits > &__os, const bitset< _Nb > &__x)
Global I/O operators for bitsets.
Definition bitset:1733
Implementation details not part of the namespace std interface.
initializer_list
A standard container for storing a fixed size sequence of elements.
Definition array:103
char_type widen(char __c) const
Widens characters.
Definition basic_ios.h:464
char_type fill() const
Retrieves the empty character.
Definition basic_ios.h:387
Template class basic_istream.
Definition istream:63
Template class basic_ostream.
Definition ostream.h:67
static constexpr bool is_integer
Definition limits:233
static constexpr int max_digits10
Definition limits:226
static constexpr int digits
Definition limits:218
static constexpr bool is_signed
Definition limits:230
static constexpr _Tp max() noexcept
Definition limits:328
static constexpr _Tp epsilon() noexcept
Definition limits:340
is_floating_point
Definition type_traits:546
common_type
Definition type_traits:2530
streamsize precision() const
Flags access.
Definition ios_base.h:765
fmtflags flags() const
Access to format flags.
Definition ios_base.h:694
A model of a linear congruential random number generator.
Definition random.h:368
static constexpr result_type multiplier
Definition random.h:383
static constexpr result_type modulus
Definition random.h:387
void seed(result_type __s=default_seed)
Reseeds the linear_congruential_engine random number generator engine sequence to the seed __s.
static constexpr result_type increment
Definition random.h:385
The Marsaglia-Zaman generator.
Definition random.h:815
void seed(result_type __sd=0u)
Seeds the initial state of the random number generator.
result_type operator()()
Gets the next random number in the sequence.
result_type operator()()
Gets the next value in the generated random number sequence.
result_type operator()()
Gets the next value in the generated random number sequence.
Produces random numbers by reordering random numbers from some base engine.
Definition random.h:1465
const _RandomNumberEngine & base() const noexcept
Definition random.h:1571
_RandomNumberEngine::result_type result_type
Definition random.h:1471
Uniform continuous distribution for random numbers.
Definition random.h:2157
param_type param() const
Returns the parameter set of the distribution.
Definition random.h:2246
A normal continuous distribution for random numbers.
Definition random.h:2394
result_type operator()(_UniformRandomNumberGenerator &__urng)
Generating functions.
Definition random.h:2513
result_type operator()(_UniformRandomNumberGenerator &__urng)
Generating functions.
Definition random.h:2975
A cauchy_distribution random number distribution.
Definition random.h:3318
param_type param() const
Returns the parameter set of the distribution.
Definition random.h:3395
result_type operator()(_UniformRandomNumberGenerator &__urng)
Generating functions.
Definition random.h:3425
A fisher_f_distribution random number distribution.
Definition random.h:3533
A Bernoulli random number distribution.
Definition random.h:4005
double p() const
Returns the p parameter of the distribution.
Definition random.h:4076
A discrete binomial random number distribution.
Definition random.h:4229
result_type operator()(_UniformRandomNumberGenerator &__urng)
Generating functions.
Definition random.h:4357
A discrete geometric random number distribution.
Definition random.h:4475
result_type operator()(_UniformRandomNumberGenerator &__urng)
Generating functions.
Definition random.h:4586
param_type param() const
Returns the parameter set of the distribution.
Definition random.h:4556
A negative_binomial_distribution random number distribution.
Definition random.h:4692
result_type operator()(_UniformRandomNumberGenerator &__urng)
Generating functions.
A discrete Poisson random number distribution.
Definition random.h:4929
result_type operator()(_UniformRandomNumberGenerator &__urng)
Generating functions.
Definition random.h:5042
friend bool operator==(const poisson_distribution &__d1, const poisson_distribution &__d2)
Return true if two Poisson distributions have the same parameters and the sequences that would be gen...
Definition random.h:5078
friend std::basic_ostream< _CharT, _Traits > & operator<<(std::basic_ostream< _CharT, _Traits > &__os, const std::poisson_distribution< _IntType1 > &__x)
Inserts a poisson_distribution random number distribution __x into the output stream __os.
friend std::basic_istream< _CharT, _Traits > & operator>>(std::basic_istream< _CharT, _Traits > &__is, std::poisson_distribution< _IntType1 > &__x)
Extracts a poisson_distribution random number distribution __x from the input stream __is.
An exponential continuous distribution for random numbers.
Definition random.h:5161
param_type param() const
Returns the parameter set of the distribution.
Definition random.h:5241
A weibull_distribution random number distribution.
Definition random.h:5383
param_type param() const
Returns the parameter set of the distribution.
Definition random.h:5463
result_type operator()(_UniformRandomNumberGenerator &__urng)
Generating functions.
Definition random.h:5493
A extreme_value_distribution random number distribution.
Definition random.h:5600
result_type operator()(_UniformRandomNumberGenerator &__urng)
Generating functions.
Definition random.h:5710
param_type param() const
Returns the parameter set of the distribution.
Definition random.h:5680
A piecewise_constant_distribution random number distribution.
Definition random.h:6070
result_type operator()(_UniformRandomNumberGenerator &__urng)
Generating functions.
Definition random.h:6242
A piecewise_linear_distribution random number distribution.
Definition random.h:6350
result_type operator()(_UniformRandomNumberGenerator &__urng)
Generating functions.
Definition random.h:6524
seed_seq() noexcept
Definition random.h:6642
uint_least32_t result_type
Definition random.h:6639
A standard container which offers fixed time access to individual elements in any order.
Definition stl_vector.h:461
constexpr iterator end() noexcept
constexpr iterator begin() noexcept
Definition stl_vector.h:988
constexpr size_type size() const noexcept
Uniform discrete distribution for random numbers. A discrete random distribution on the range with e...
param_type param() const
Returns the parameter set of the distribution.