libstdc++
ext/random.tcc
Go to the documentation of this file.
1 // Random number extensions -*- C++ -*-
2 
3 // Copyright (C) 2012-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 ext/random.tcc
26  * This is an internal header file, included by other library headers.
27  * Do not attempt to use it directly. @headername{ext/random}
28  */
29 
30 #ifndef _EXT_RANDOM_TCC
31 #define _EXT_RANDOM_TCC 1
32 
33 #ifdef _GLIBCXX_SYSHDR
34 #pragma GCC system_header
35 #endif
36 
37 #include <bits/requires_hosted.h> // GNU extensions are currently omitted
38 
39 namespace __gnu_cxx _GLIBCXX_VISIBILITY(default)
40 {
41 _GLIBCXX_BEGIN_NAMESPACE_VERSION
42 
43 #if __BYTE_ORDER__ == __ORDER_LITTLE_ENDIAN__
44 
45  template<typename _UIntType, size_t __m,
46  size_t __pos1, size_t __sl1, size_t __sl2,
47  size_t __sr1, size_t __sr2,
48  uint32_t __msk1, uint32_t __msk2,
49  uint32_t __msk3, uint32_t __msk4,
50  uint32_t __parity1, uint32_t __parity2,
51  uint32_t __parity3, uint32_t __parity4>
52  void simd_fast_mersenne_twister_engine<_UIntType, __m,
53  __pos1, __sl1, __sl2, __sr1, __sr2,
54  __msk1, __msk2, __msk3, __msk4,
55  __parity1, __parity2, __parity3,
56  __parity4>::
57  seed(_UIntType __seed)
58  {
59  _M_state32[0] = static_cast<uint32_t>(__seed);
60  for (size_t __i = 1; __i < _M_nstate32; ++__i)
61  _M_state32[__i] = (1812433253UL
62  * (_M_state32[__i - 1] ^ (_M_state32[__i - 1] >> 30))
63  + __i);
64  _M_pos = state_size;
65  _M_period_certification();
66  }
67 
68 
69  namespace {
70 
71  inline uint32_t _Func1(uint32_t __x)
72  {
73  return (__x ^ (__x >> 27)) * UINT32_C(1664525);
74  }
75 
76  inline uint32_t _Func2(uint32_t __x)
77  {
78  return (__x ^ (__x >> 27)) * UINT32_C(1566083941);
79  }
80 
81  }
82 
83 
84  template<typename _UIntType, size_t __m,
85  size_t __pos1, size_t __sl1, size_t __sl2,
86  size_t __sr1, size_t __sr2,
87  uint32_t __msk1, uint32_t __msk2,
88  uint32_t __msk3, uint32_t __msk4,
89  uint32_t __parity1, uint32_t __parity2,
90  uint32_t __parity3, uint32_t __parity4>
91  template<typename _Sseq>
92  auto
93  simd_fast_mersenne_twister_engine<_UIntType, __m,
94  __pos1, __sl1, __sl2, __sr1, __sr2,
95  __msk1, __msk2, __msk3, __msk4,
96  __parity1, __parity2, __parity3,
97  __parity4>::
98  seed(_Sseq& __q)
99  -> _If_seed_seq<_Sseq>
100  {
101  size_t __lag;
102 
103  if (_M_nstate32 >= 623)
104  __lag = 11;
105  else if (_M_nstate32 >= 68)
106  __lag = 7;
107  else if (_M_nstate32 >= 39)
108  __lag = 5;
109  else
110  __lag = 3;
111  const size_t __mid = (_M_nstate32 - __lag) / 2;
112 
113  std::fill(_M_state32, _M_state32 + _M_nstate32, UINT32_C(0x8b8b8b8b));
114  uint32_t __arr[_M_nstate32];
115  __q.generate(__arr + 0, __arr + _M_nstate32);
116 
117  uint32_t __r = _Func1(_M_state32[0] ^ _M_state32[__mid]
118  ^ _M_state32[_M_nstate32 - 1]);
119  _M_state32[__mid] += __r;
120  __r += _M_nstate32;
121  _M_state32[__mid + __lag] += __r;
122  _M_state32[0] = __r;
123 
124  for (size_t __i = 1, __j = 0; __j < _M_nstate32; ++__j)
125  {
126  __r = _Func1(_M_state32[__i]
127  ^ _M_state32[(__i + __mid) % _M_nstate32]
128  ^ _M_state32[(__i + _M_nstate32 - 1) % _M_nstate32]);
129  _M_state32[(__i + __mid) % _M_nstate32] += __r;
130  __r += __arr[__j] + __i;
131  _M_state32[(__i + __mid + __lag) % _M_nstate32] += __r;
132  _M_state32[__i] = __r;
133  __i = (__i + 1) % _M_nstate32;
134  }
135  for (size_t __j = 0; __j < _M_nstate32; ++__j)
136  {
137  const size_t __i = (__j + 1) % _M_nstate32;
138  __r = _Func2(_M_state32[__i]
139  + _M_state32[(__i + __mid) % _M_nstate32]
140  + _M_state32[(__i + _M_nstate32 - 1) % _M_nstate32]);
141  _M_state32[(__i + __mid) % _M_nstate32] ^= __r;
142  __r -= __i;
143  _M_state32[(__i + __mid + __lag) % _M_nstate32] ^= __r;
144  _M_state32[__i] = __r;
145  }
146 
147  _M_pos = state_size;
148  _M_period_certification();
149  }
150 
151 
152  template<typename _UIntType, size_t __m,
153  size_t __pos1, size_t __sl1, size_t __sl2,
154  size_t __sr1, size_t __sr2,
155  uint32_t __msk1, uint32_t __msk2,
156  uint32_t __msk3, uint32_t __msk4,
157  uint32_t __parity1, uint32_t __parity2,
158  uint32_t __parity3, uint32_t __parity4>
159  void simd_fast_mersenne_twister_engine<_UIntType, __m,
160  __pos1, __sl1, __sl2, __sr1, __sr2,
161  __msk1, __msk2, __msk3, __msk4,
162  __parity1, __parity2, __parity3,
163  __parity4>::
164  _M_period_certification(void)
165  {
166  static const uint32_t __parity[4] = { __parity1, __parity2,
167  __parity3, __parity4 };
168  uint32_t __inner = 0;
169  for (size_t __i = 0; __i < 4; ++__i)
170  if (__parity[__i] != 0)
171  __inner ^= _M_state32[__i] & __parity[__i];
172 
173  if (__builtin_parity(__inner) & 1)
174  return;
175  for (size_t __i = 0; __i < 4; ++__i)
176  if (__parity[__i] != 0)
177  {
178  _M_state32[__i] ^= 1 << (__builtin_ffs(__parity[__i]) - 1);
179  return;
180  }
181  __builtin_unreachable();
182  }
183 
184 
185  template<typename _UIntType, size_t __m,
186  size_t __pos1, size_t __sl1, size_t __sl2,
187  size_t __sr1, size_t __sr2,
188  uint32_t __msk1, uint32_t __msk2,
189  uint32_t __msk3, uint32_t __msk4,
190  uint32_t __parity1, uint32_t __parity2,
191  uint32_t __parity3, uint32_t __parity4>
192  void simd_fast_mersenne_twister_engine<_UIntType, __m,
193  __pos1, __sl1, __sl2, __sr1, __sr2,
194  __msk1, __msk2, __msk3, __msk4,
195  __parity1, __parity2, __parity3,
196  __parity4>::
197  discard(unsigned long long __z)
198  {
199  while (__z > state_size - _M_pos)
200  {
201  __z -= state_size - _M_pos;
202 
203  _M_gen_rand();
204  }
205 
206  _M_pos += __z;
207  }
208 
209 
210 #ifndef _GLIBCXX_OPT_HAVE_RANDOM_SFMT_GEN_READ
211 
212  namespace {
213 
214  template<size_t __shift>
215  inline void __rshift(uint32_t *__out, const uint32_t *__in)
216  {
217  uint64_t __th = ((static_cast<uint64_t>(__in[3]) << 32)
218  | static_cast<uint64_t>(__in[2]));
219  uint64_t __tl = ((static_cast<uint64_t>(__in[1]) << 32)
220  | static_cast<uint64_t>(__in[0]));
221 
222  uint64_t __oh = __th >> (__shift * 8);
223  uint64_t __ol = __tl >> (__shift * 8);
224  __ol |= __th << (64 - __shift * 8);
225  __out[1] = static_cast<uint32_t>(__ol >> 32);
226  __out[0] = static_cast<uint32_t>(__ol);
227  __out[3] = static_cast<uint32_t>(__oh >> 32);
228  __out[2] = static_cast<uint32_t>(__oh);
229  }
230 
231 
232  template<size_t __shift>
233  inline void __lshift(uint32_t *__out, const uint32_t *__in)
234  {
235  uint64_t __th = ((static_cast<uint64_t>(__in[3]) << 32)
236  | static_cast<uint64_t>(__in[2]));
237  uint64_t __tl = ((static_cast<uint64_t>(__in[1]) << 32)
238  | static_cast<uint64_t>(__in[0]));
239 
240  uint64_t __oh = __th << (__shift * 8);
241  uint64_t __ol = __tl << (__shift * 8);
242  __oh |= __tl >> (64 - __shift * 8);
243  __out[1] = static_cast<uint32_t>(__ol >> 32);
244  __out[0] = static_cast<uint32_t>(__ol);
245  __out[3] = static_cast<uint32_t>(__oh >> 32);
246  __out[2] = static_cast<uint32_t>(__oh);
247  }
248 
249 
250  template<size_t __sl1, size_t __sl2, size_t __sr1, size_t __sr2,
251  uint32_t __msk1, uint32_t __msk2, uint32_t __msk3, uint32_t __msk4>
252  inline void __recursion(uint32_t *__r,
253  const uint32_t *__a, const uint32_t *__b,
254  const uint32_t *__c, const uint32_t *__d)
255  {
256  uint32_t __x[4];
257  uint32_t __y[4];
258 
259  __lshift<__sl2>(__x, __a);
260  __rshift<__sr2>(__y, __c);
261  __r[0] = (__a[0] ^ __x[0] ^ ((__b[0] >> __sr1) & __msk1)
262  ^ __y[0] ^ (__d[0] << __sl1));
263  __r[1] = (__a[1] ^ __x[1] ^ ((__b[1] >> __sr1) & __msk2)
264  ^ __y[1] ^ (__d[1] << __sl1));
265  __r[2] = (__a[2] ^ __x[2] ^ ((__b[2] >> __sr1) & __msk3)
266  ^ __y[2] ^ (__d[2] << __sl1));
267  __r[3] = (__a[3] ^ __x[3] ^ ((__b[3] >> __sr1) & __msk4)
268  ^ __y[3] ^ (__d[3] << __sl1));
269  }
270 
271  }
272 
273 
274  template<typename _UIntType, size_t __m,
275  size_t __pos1, size_t __sl1, size_t __sl2,
276  size_t __sr1, size_t __sr2,
277  uint32_t __msk1, uint32_t __msk2,
278  uint32_t __msk3, uint32_t __msk4,
279  uint32_t __parity1, uint32_t __parity2,
280  uint32_t __parity3, uint32_t __parity4>
281  void simd_fast_mersenne_twister_engine<_UIntType, __m,
282  __pos1, __sl1, __sl2, __sr1, __sr2,
283  __msk1, __msk2, __msk3, __msk4,
284  __parity1, __parity2, __parity3,
285  __parity4>::
286  _M_gen_rand(void)
287  {
288  const uint32_t *__r1 = &_M_state32[_M_nstate32 - 8];
289  const uint32_t *__r2 = &_M_state32[_M_nstate32 - 4];
290  static constexpr size_t __pos1_32 = __pos1 * 4;
291 
292  size_t __i;
293  for (__i = 0; __i < _M_nstate32 - __pos1_32; __i += 4)
294  {
295  __recursion<__sl1, __sl2, __sr1, __sr2,
296  __msk1, __msk2, __msk3, __msk4>
297  (&_M_state32[__i], &_M_state32[__i],
298  &_M_state32[__i + __pos1_32], __r1, __r2);
299  __r1 = __r2;
300  __r2 = &_M_state32[__i];
301  }
302 
303  for (; __i < _M_nstate32; __i += 4)
304  {
305  __recursion<__sl1, __sl2, __sr1, __sr2,
306  __msk1, __msk2, __msk3, __msk4>
307  (&_M_state32[__i], &_M_state32[__i],
308  &_M_state32[__i + __pos1_32 - _M_nstate32], __r1, __r2);
309  __r1 = __r2;
310  __r2 = &_M_state32[__i];
311  }
312 
313  _M_pos = 0;
314  }
315 
316 #endif
317 
318 #ifndef _GLIBCXX_OPT_HAVE_RANDOM_SFMT_OPERATOREQUAL
319  template<typename _UIntType, size_t __m,
320  size_t __pos1, size_t __sl1, size_t __sl2,
321  size_t __sr1, size_t __sr2,
322  uint32_t __msk1, uint32_t __msk2,
323  uint32_t __msk3, uint32_t __msk4,
324  uint32_t __parity1, uint32_t __parity2,
325  uint32_t __parity3, uint32_t __parity4>
326  bool
327  operator==(const __gnu_cxx::simd_fast_mersenne_twister_engine<_UIntType,
328  __m, __pos1, __sl1, __sl2, __sr1, __sr2,
329  __msk1, __msk2, __msk3, __msk4,
330  __parity1, __parity2, __parity3, __parity4>& __lhs,
331  const __gnu_cxx::simd_fast_mersenne_twister_engine<_UIntType,
332  __m, __pos1, __sl1, __sl2, __sr1, __sr2,
333  __msk1, __msk2, __msk3, __msk4,
334  __parity1, __parity2, __parity3, __parity4>& __rhs)
335  {
336  typedef __gnu_cxx::simd_fast_mersenne_twister_engine<_UIntType,
337  __m, __pos1, __sl1, __sl2, __sr1, __sr2,
338  __msk1, __msk2, __msk3, __msk4,
339  __parity1, __parity2, __parity3, __parity4> __engine;
340  return (std::equal(__lhs._M_stateT,
341  __lhs._M_stateT + __engine::state_size,
342  __rhs._M_stateT)
343  && __lhs._M_pos == __rhs._M_pos);
344  }
345 #endif
346 
347  template<typename _UIntType, size_t __m,
348  size_t __pos1, size_t __sl1, size_t __sl2,
349  size_t __sr1, size_t __sr2,
350  uint32_t __msk1, uint32_t __msk2,
351  uint32_t __msk3, uint32_t __msk4,
352  uint32_t __parity1, uint32_t __parity2,
353  uint32_t __parity3, uint32_t __parity4,
354  typename _CharT, typename _Traits>
357  const __gnu_cxx::simd_fast_mersenne_twister_engine<_UIntType,
358  __m, __pos1, __sl1, __sl2, __sr1, __sr2,
359  __msk1, __msk2, __msk3, __msk4,
360  __parity1, __parity2, __parity3, __parity4>& __x)
361  {
362  typedef std::basic_ostream<_CharT, _Traits> __ostream_type;
363  typedef typename __ostream_type::ios_base __ios_base;
364 
365  const typename __ios_base::fmtflags __flags = __os.flags();
366  const _CharT __fill = __os.fill();
367  const _CharT __space = __os.widen(' ');
369  __os.fill(__space);
370 
371  for (size_t __i = 0; __i < __x._M_nstate32; ++__i)
372  __os << __x._M_state32[__i] << __space;
373  __os << __x._M_pos;
374 
375  __os.flags(__flags);
376  __os.fill(__fill);
377  return __os;
378  }
379 
380 
381  template<typename _UIntType, size_t __m,
382  size_t __pos1, size_t __sl1, size_t __sl2,
383  size_t __sr1, size_t __sr2,
384  uint32_t __msk1, uint32_t __msk2,
385  uint32_t __msk3, uint32_t __msk4,
386  uint32_t __parity1, uint32_t __parity2,
387  uint32_t __parity3, uint32_t __parity4,
388  typename _CharT, typename _Traits>
391  __gnu_cxx::simd_fast_mersenne_twister_engine<_UIntType,
392  __m, __pos1, __sl1, __sl2, __sr1, __sr2,
393  __msk1, __msk2, __msk3, __msk4,
394  __parity1, __parity2, __parity3, __parity4>& __x)
395  {
396  typedef std::basic_istream<_CharT, _Traits> __istream_type;
397  typedef typename __istream_type::ios_base __ios_base;
398 
399  const typename __ios_base::fmtflags __flags = __is.flags();
401 
402  for (size_t __i = 0; __i < __x._M_nstate32; ++__i)
403  __is >> __x._M_state32[__i];
404  __is >> __x._M_pos;
405 
406  __is.flags(__flags);
407  return __is;
408  }
409 
410 #endif // __BYTE_ORDER__ == __ORDER_LITTLE_ENDIAN__
411 
412  /**
413  * Iteration method due to M.D. J<o:>hnk.
414  *
415  * M.D. J<o:>hnk, Erzeugung von betaverteilten und gammaverteilten
416  * Zufallszahlen, Metrika, Volume 8, 1964
417  */
418  template<typename _RealType>
419  template<typename _UniformRandomNumberGenerator>
420  typename beta_distribution<_RealType>::result_type
421  beta_distribution<_RealType>::
422  operator()(_UniformRandomNumberGenerator& __urng,
423  const param_type& __param)
424  {
425  std::__detail::_Adaptor<_UniformRandomNumberGenerator, result_type>
426  __aurng(__urng);
427 
428  result_type __x, __y;
429  do
430  {
431  __x = std::exp(std::log(__aurng()) / __param.alpha());
432  __y = std::exp(std::log(__aurng()) / __param.beta());
433  }
434  while (__x + __y > result_type(1));
435 
436  return __x / (__x + __y);
437  }
438 
439  template<typename _RealType>
440  template<typename _OutputIterator,
441  typename _UniformRandomNumberGenerator>
442  void
443  beta_distribution<_RealType>::
444  __generate_impl(_OutputIterator __f, _OutputIterator __t,
445  _UniformRandomNumberGenerator& __urng,
446  const param_type& __param)
447  {
448  __glibcxx_function_requires(_OutputIteratorConcept<_OutputIterator,
449  result_type>)
450 
451  std::__detail::_Adaptor<_UniformRandomNumberGenerator, result_type>
452  __aurng(__urng);
453 
454  while (__f != __t)
455  {
456  result_type __x, __y;
457  do
458  {
459  __x = std::exp(std::log(__aurng()) / __param.alpha());
460  __y = std::exp(std::log(__aurng()) / __param.beta());
461  }
462  while (__x + __y > result_type(1));
463 
464  *__f++ = __x / (__x + __y);
465  }
466  }
467 
468  template<typename _RealType, typename _CharT, typename _Traits>
471  const __gnu_cxx::beta_distribution<_RealType>& __x)
472  {
473  typedef std::basic_ostream<_CharT, _Traits> __ostream_type;
474  typedef typename __ostream_type::ios_base __ios_base;
475 
476  const typename __ios_base::fmtflags __flags = __os.flags();
477  const _CharT __fill = __os.fill();
478  const std::streamsize __precision = __os.precision();
479  const _CharT __space = __os.widen(' ');
481  __os.fill(__space);
483 
484  __os << __x.alpha() << __space << __x.beta();
485 
486  __os.flags(__flags);
487  __os.fill(__fill);
488  __os.precision(__precision);
489  return __os;
490  }
491 
492  template<typename _RealType, typename _CharT, typename _Traits>
495  __gnu_cxx::beta_distribution<_RealType>& __x)
496  {
497  typedef std::basic_istream<_CharT, _Traits> __istream_type;
498  typedef typename __istream_type::ios_base __ios_base;
499 
500  const typename __ios_base::fmtflags __flags = __is.flags();
502 
503  _RealType __alpha_val, __beta_val;
504  __is >> __alpha_val >> __beta_val;
505  __x.param(typename __gnu_cxx::beta_distribution<_RealType>::
506  param_type(__alpha_val, __beta_val));
507 
508  __is.flags(__flags);
509  return __is;
510  }
511 
512 
513  template<std::size_t _Dimen, typename _RealType>
514  template<typename _InputIterator1, typename _InputIterator2>
515  void
516  normal_mv_distribution<_Dimen, _RealType>::param_type::
517  _M_init_full(_InputIterator1 __meanbegin, _InputIterator1 __meanend,
518  _InputIterator2 __varcovbegin, _InputIterator2 __varcovend)
519  {
520  __glibcxx_function_requires(_InputIteratorConcept<_InputIterator1>)
521  __glibcxx_function_requires(_InputIteratorConcept<_InputIterator2>)
522  std::fill(std::copy(__meanbegin, __meanend, _M_mean.begin()),
523  _M_mean.end(), _RealType(0));
524 
525  // Perform the Cholesky decomposition
526  auto __w = _M_t.begin();
527  for (size_t __j = 0; __j < _Dimen; ++__j)
528  {
529  _RealType __sum = _RealType(0);
530 
531  auto __slitbegin = __w;
532  auto __cit = _M_t.begin();
533  for (size_t __i = 0; __i < __j; ++__i)
534  {
535  auto __slit = __slitbegin;
536  _RealType __s = *__varcovbegin++;
537  for (size_t __k = 0; __k < __i; ++__k)
538  __s -= *__slit++ * *__cit++;
539 
540  *__w++ = __s /= *__cit++;
541  __sum += __s * __s;
542  }
543 
544  __sum = *__varcovbegin - __sum;
545  if (__builtin_expect(__sum <= _RealType(0), 0))
546  std::__throw_runtime_error(__N("normal_mv_distribution::"
547  "param_type::_M_init_full"));
548  *__w++ = std::sqrt(__sum);
549 
550  std::advance(__varcovbegin, _Dimen - __j);
551  }
552  }
553 
554  template<std::size_t _Dimen, typename _RealType>
555  template<typename _InputIterator1, typename _InputIterator2>
556  void
557  normal_mv_distribution<_Dimen, _RealType>::param_type::
558  _M_init_lower(_InputIterator1 __meanbegin, _InputIterator1 __meanend,
559  _InputIterator2 __varcovbegin, _InputIterator2 __varcovend)
560  {
561  __glibcxx_function_requires(_InputIteratorConcept<_InputIterator1>)
562  __glibcxx_function_requires(_InputIteratorConcept<_InputIterator2>)
563  std::fill(std::copy(__meanbegin, __meanend, _M_mean.begin()),
564  _M_mean.end(), _RealType(0));
565 
566  // Perform the Cholesky decomposition
567  auto __w = _M_t.begin();
568  for (size_t __j = 0; __j < _Dimen; ++__j)
569  {
570  _RealType __sum = _RealType(0);
571 
572  auto __slitbegin = __w;
573  auto __cit = _M_t.begin();
574  for (size_t __i = 0; __i < __j; ++__i)
575  {
576  auto __slit = __slitbegin;
577  _RealType __s = *__varcovbegin++;
578  for (size_t __k = 0; __k < __i; ++__k)
579  __s -= *__slit++ * *__cit++;
580 
581  *__w++ = __s /= *__cit++;
582  __sum += __s * __s;
583  }
584 
585  __sum = *__varcovbegin++ - __sum;
586  if (__builtin_expect(__sum <= _RealType(0), 0))
587  std::__throw_runtime_error(__N("normal_mv_distribution::"
588  "param_type::_M_init_lower"));
589  *__w++ = std::sqrt(__sum);
590  }
591  }
592 
593  template<std::size_t _Dimen, typename _RealType>
594  template<typename _InputIterator1, typename _InputIterator2>
595  void
596  normal_mv_distribution<_Dimen, _RealType>::param_type::
597  _M_init_diagonal(_InputIterator1 __meanbegin, _InputIterator1 __meanend,
598  _InputIterator2 __varbegin, _InputIterator2 __varend)
599  {
600  __glibcxx_function_requires(_InputIteratorConcept<_InputIterator1>)
601  __glibcxx_function_requires(_InputIteratorConcept<_InputIterator2>)
602  std::fill(std::copy(__meanbegin, __meanend, _M_mean.begin()),
603  _M_mean.end(), _RealType(0));
604 
605  auto __w = _M_t.begin();
606  size_t __step = 0;
607  while (__varbegin != __varend)
608  {
609  std::fill_n(__w, __step, _RealType(0));
610  __w += __step++;
611  if (__builtin_expect(*__varbegin < _RealType(0), 0))
612  std::__throw_runtime_error(__N("normal_mv_distribution::"
613  "param_type::_M_init_diagonal"));
614  *__w++ = std::sqrt(*__varbegin++);
615  }
616  }
617 
618  template<std::size_t _Dimen, typename _RealType>
619  template<typename _UniformRandomNumberGenerator>
620  typename normal_mv_distribution<_Dimen, _RealType>::result_type
621  normal_mv_distribution<_Dimen, _RealType>::
622  operator()(_UniformRandomNumberGenerator& __urng,
623  const param_type& __param)
624  {
625  result_type __ret;
626 
627  _M_nd.__generate(__ret.begin(), __ret.end(), __urng);
628 
629  auto __t_it = __param._M_t.crbegin();
630  for (size_t __i = _Dimen; __i > 0; --__i)
631  {
632  _RealType __sum = _RealType(0);
633  for (size_t __j = __i; __j > 0; --__j)
634  __sum += __ret[__j - 1] * *__t_it++;
635  __ret[__i - 1] = __sum;
636  }
637 
638  return __ret;
639  }
640 
641  template<std::size_t _Dimen, typename _RealType>
642  template<typename _ForwardIterator, typename _UniformRandomNumberGenerator>
643  void
644  normal_mv_distribution<_Dimen, _RealType>::
645  __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
646  _UniformRandomNumberGenerator& __urng,
647  const param_type& __param)
648  {
649  __glibcxx_function_requires(_Mutable_ForwardIteratorConcept<
650  _ForwardIterator>)
651  while (__f != __t)
652  *__f++ = this->operator()(__urng, __param);
653  }
654 
655  template<size_t _Dimen, typename _RealType>
656  bool
657  operator==(const __gnu_cxx::normal_mv_distribution<_Dimen, _RealType>&
658  __d1,
659  const __gnu_cxx::normal_mv_distribution<_Dimen, _RealType>&
660  __d2)
661  {
662  return __d1._M_param == __d2._M_param && __d1._M_nd == __d2._M_nd;
663  }
664 
665  template<size_t _Dimen, typename _RealType, typename _CharT, typename _Traits>
668  const __gnu_cxx::normal_mv_distribution<_Dimen, _RealType>& __x)
669  {
670  typedef std::basic_ostream<_CharT, _Traits> __ostream_type;
671  typedef typename __ostream_type::ios_base __ios_base;
672 
673  const typename __ios_base::fmtflags __flags = __os.flags();
674  const _CharT __fill = __os.fill();
675  const std::streamsize __precision = __os.precision();
676  const _CharT __space = __os.widen(' ');
678  __os.fill(__space);
680 
681  auto __mean = __x._M_param.mean();
682  for (auto __it : __mean)
683  __os << __it << __space;
684  auto __t = __x._M_param.varcov();
685  for (auto __it : __t)
686  __os << __it << __space;
687 
688  __os << __x._M_nd;
689 
690  __os.flags(__flags);
691  __os.fill(__fill);
692  __os.precision(__precision);
693  return __os;
694  }
695 
696  template<size_t _Dimen, typename _RealType, typename _CharT, typename _Traits>
699  __gnu_cxx::normal_mv_distribution<_Dimen, _RealType>& __x)
700  {
701  typedef std::basic_istream<_CharT, _Traits> __istream_type;
702  typedef typename __istream_type::ios_base __ios_base;
703 
704  const typename __ios_base::fmtflags __flags = __is.flags();
706 
708  for (auto& __it : __mean)
709  __is >> __it;
710  std::array<_RealType, _Dimen * (_Dimen + 1) / 2> __varcov;
711  for (auto& __it : __varcov)
712  __is >> __it;
713 
714  __is >> __x._M_nd;
715 
716  // The param_type temporary is built with a private constructor,
717  // to skip the Cholesky decomposition that would be performed
718  // otherwise.
719  __x.param(typename normal_mv_distribution<_Dimen, _RealType>::
720  param_type(__mean, __varcov));
721 
722  __is.flags(__flags);
723  return __is;
724  }
725 
726 
727  template<typename _RealType>
728  template<typename _OutputIterator,
729  typename _UniformRandomNumberGenerator>
730  void
731  rice_distribution<_RealType>::
732  __generate_impl(_OutputIterator __f, _OutputIterator __t,
733  _UniformRandomNumberGenerator& __urng,
734  const param_type& __p)
735  {
736  __glibcxx_function_requires(_OutputIteratorConcept<_OutputIterator,
737  result_type>)
738 
739  while (__f != __t)
740  {
742  __px(__p.nu(), __p.sigma()), __py(result_type(0), __p.sigma());
743  result_type __x = this->_M_ndx(__px, __urng);
744  result_type __y = this->_M_ndy(__py, __urng);
745 #if _GLIBCXX_USE_C99_MATH_FUNCS
746  *__f++ = std::hypot(__x, __y);
747 #else
748  *__f++ = std::sqrt(__x * __x + __y * __y);
749 #endif
750  }
751  }
752 
753  template<typename _RealType, typename _CharT, typename _Traits>
756  const rice_distribution<_RealType>& __x)
757  {
758  typedef std::basic_ostream<_CharT, _Traits> __ostream_type;
759  typedef typename __ostream_type::ios_base __ios_base;
760 
761  const typename __ios_base::fmtflags __flags = __os.flags();
762  const _CharT __fill = __os.fill();
763  const std::streamsize __precision = __os.precision();
764  const _CharT __space = __os.widen(' ');
766  __os.fill(__space);
768 
769  __os << __x.nu() << __space << __x.sigma();
770  __os << __space << __x._M_ndx;
771  __os << __space << __x._M_ndy;
772 
773  __os.flags(__flags);
774  __os.fill(__fill);
775  __os.precision(__precision);
776  return __os;
777  }
778 
779  template<typename _RealType, typename _CharT, typename _Traits>
782  rice_distribution<_RealType>& __x)
783  {
784  typedef std::basic_istream<_CharT, _Traits> __istream_type;
785  typedef typename __istream_type::ios_base __ios_base;
786 
787  const typename __ios_base::fmtflags __flags = __is.flags();
789 
790  _RealType __nu_val, __sigma_val;
791  __is >> __nu_val >> __sigma_val;
792  __is >> __x._M_ndx;
793  __is >> __x._M_ndy;
794  __x.param(typename rice_distribution<_RealType>::
795  param_type(__nu_val, __sigma_val));
796 
797  __is.flags(__flags);
798  return __is;
799  }
800 
801 
802  template<typename _RealType>
803  template<typename _OutputIterator,
804  typename _UniformRandomNumberGenerator>
805  void
806  nakagami_distribution<_RealType>::
807  __generate_impl(_OutputIterator __f, _OutputIterator __t,
808  _UniformRandomNumberGenerator& __urng,
809  const param_type& __p)
810  {
811  __glibcxx_function_requires(_OutputIteratorConcept<_OutputIterator,
812  result_type>)
813 
814  typename std::gamma_distribution<result_type>::param_type
815  __pg(__p.mu(), __p.omega() / __p.mu());
816  while (__f != __t)
817  *__f++ = std::sqrt(this->_M_gd(__pg, __urng));
818  }
819 
820  template<typename _RealType, typename _CharT, typename _Traits>
821  std::basic_ostream<_CharT, _Traits>&
822  operator<<(std::basic_ostream<_CharT, _Traits>& __os,
823  const nakagami_distribution<_RealType>& __x)
824  {
825  typedef std::basic_ostream<_CharT, _Traits> __ostream_type;
826  typedef typename __ostream_type::ios_base __ios_base;
827 
828  const typename __ios_base::fmtflags __flags = __os.flags();
829  const _CharT __fill = __os.fill();
830  const std::streamsize __precision = __os.precision();
831  const _CharT __space = __os.widen(' ');
833  __os.fill(__space);
835 
836  __os << __x.mu() << __space << __x.omega();
837  __os << __space << __x._M_gd;
838 
839  __os.flags(__flags);
840  __os.fill(__fill);
841  __os.precision(__precision);
842  return __os;
843  }
844 
845  template<typename _RealType, typename _CharT, typename _Traits>
848  nakagami_distribution<_RealType>& __x)
849  {
850  typedef std::basic_istream<_CharT, _Traits> __istream_type;
851  typedef typename __istream_type::ios_base __ios_base;
852 
853  const typename __ios_base::fmtflags __flags = __is.flags();
855 
856  _RealType __mu_val, __omega_val;
857  __is >> __mu_val >> __omega_val;
858  __is >> __x._M_gd;
859  __x.param(typename nakagami_distribution<_RealType>::
860  param_type(__mu_val, __omega_val));
861 
862  __is.flags(__flags);
863  return __is;
864  }
865 
866 
867  template<typename _RealType>
868  template<typename _OutputIterator,
869  typename _UniformRandomNumberGenerator>
870  void
871  pareto_distribution<_RealType>::
872  __generate_impl(_OutputIterator __f, _OutputIterator __t,
873  _UniformRandomNumberGenerator& __urng,
874  const param_type& __p)
875  {
876  __glibcxx_function_requires(_OutputIteratorConcept<_OutputIterator,
877  result_type>)
878 
879  result_type __mu_val = __p.mu();
880  result_type __malphinv = -result_type(1) / __p.alpha();
881  while (__f != __t)
882  *__f++ = __mu_val * std::pow(this->_M_ud(__urng), __malphinv);
883  }
884 
885  template<typename _RealType, typename _CharT, typename _Traits>
886  std::basic_ostream<_CharT, _Traits>&
887  operator<<(std::basic_ostream<_CharT, _Traits>& __os,
888  const pareto_distribution<_RealType>& __x)
889  {
890  typedef std::basic_ostream<_CharT, _Traits> __ostream_type;
891  typedef typename __ostream_type::ios_base __ios_base;
892 
893  const typename __ios_base::fmtflags __flags = __os.flags();
894  const _CharT __fill = __os.fill();
895  const std::streamsize __precision = __os.precision();
896  const _CharT __space = __os.widen(' ');
898  __os.fill(__space);
900 
901  __os << __x.alpha() << __space << __x.mu();
902  __os << __space << __x._M_ud;
903 
904  __os.flags(__flags);
905  __os.fill(__fill);
906  __os.precision(__precision);
907  return __os;
908  }
909 
910  template<typename _RealType, typename _CharT, typename _Traits>
913  pareto_distribution<_RealType>& __x)
914  {
915  typedef std::basic_istream<_CharT, _Traits> __istream_type;
916  typedef typename __istream_type::ios_base __ios_base;
917 
918  const typename __ios_base::fmtflags __flags = __is.flags();
920 
921  _RealType __alpha_val, __mu_val;
922  __is >> __alpha_val >> __mu_val;
923  __is >> __x._M_ud;
924  __x.param(typename pareto_distribution<_RealType>::
925  param_type(__alpha_val, __mu_val));
926 
927  __is.flags(__flags);
928  return __is;
929  }
930 
931 
932  template<typename _RealType>
933  template<typename _UniformRandomNumberGenerator>
934  typename k_distribution<_RealType>::result_type
935  k_distribution<_RealType>::
936  operator()(_UniformRandomNumberGenerator& __urng)
937  {
938  result_type __x = this->_M_gd1(__urng);
939  result_type __y = this->_M_gd2(__urng);
940  return std::sqrt(__x * __y);
941  }
942 
943  template<typename _RealType>
944  template<typename _UniformRandomNumberGenerator>
945  typename k_distribution<_RealType>::result_type
946  k_distribution<_RealType>::
947  operator()(_UniformRandomNumberGenerator& __urng,
948  const param_type& __p)
949  {
951  __p1(__p.lambda(), result_type(1) / __p.lambda()),
952  __p2(__p.nu(), __p.mu() / __p.nu());
953  result_type __x = this->_M_gd1(__p1, __urng);
954  result_type __y = this->_M_gd2(__p2, __urng);
955  return std::sqrt(__x * __y);
956  }
957 
958  template<typename _RealType>
959  template<typename _OutputIterator,
960  typename _UniformRandomNumberGenerator>
961  void
962  k_distribution<_RealType>::
963  __generate_impl(_OutputIterator __f, _OutputIterator __t,
964  _UniformRandomNumberGenerator& __urng,
965  const param_type& __p)
966  {
967  __glibcxx_function_requires(_OutputIteratorConcept<_OutputIterator,
968  result_type>)
969 
970  typename std::gamma_distribution<result_type>::param_type
971  __p1(__p.lambda(), result_type(1) / __p.lambda()),
972  __p2(__p.nu(), __p.mu() / __p.nu());
973  while (__f != __t)
974  {
975  result_type __x = this->_M_gd1(__p1, __urng);
976  result_type __y = this->_M_gd2(__p2, __urng);
977  *__f++ = std::sqrt(__x * __y);
978  }
979  }
980 
981  template<typename _RealType, typename _CharT, typename _Traits>
984  const k_distribution<_RealType>& __x)
985  {
986  typedef std::basic_ostream<_CharT, _Traits> __ostream_type;
987  typedef typename __ostream_type::ios_base __ios_base;
988 
989  const typename __ios_base::fmtflags __flags = __os.flags();
990  const _CharT __fill = __os.fill();
991  const std::streamsize __precision = __os.precision();
992  const _CharT __space = __os.widen(' ');
994  __os.fill(__space);
996 
997  __os << __x.lambda() << __space << __x.mu() << __space << __x.nu();
998  __os << __space << __x._M_gd1;
999  __os << __space << __x._M_gd2;
1000 
1001  __os.flags(__flags);
1002  __os.fill(__fill);
1003  __os.precision(__precision);
1004  return __os;
1005  }
1006 
1007  template<typename _RealType, typename _CharT, typename _Traits>
1010  k_distribution<_RealType>& __x)
1011  {
1012  typedef std::basic_istream<_CharT, _Traits> __istream_type;
1013  typedef typename __istream_type::ios_base __ios_base;
1014 
1015  const typename __ios_base::fmtflags __flags = __is.flags();
1017 
1018  _RealType __lambda_val, __mu_val, __nu_val;
1019  __is >> __lambda_val >> __mu_val >> __nu_val;
1020  __is >> __x._M_gd1;
1021  __is >> __x._M_gd2;
1022  __x.param(typename k_distribution<_RealType>::
1023  param_type(__lambda_val, __mu_val, __nu_val));
1024 
1025  __is.flags(__flags);
1026  return __is;
1027  }
1028 
1029 
1030  template<typename _RealType>
1031  template<typename _OutputIterator,
1032  typename _UniformRandomNumberGenerator>
1033  void
1034  arcsine_distribution<_RealType>::
1035  __generate_impl(_OutputIterator __f, _OutputIterator __t,
1036  _UniformRandomNumberGenerator& __urng,
1037  const param_type& __p)
1038  {
1039  __glibcxx_function_requires(_OutputIteratorConcept<_OutputIterator,
1040  result_type>)
1041 
1042  result_type __dif = __p.b() - __p.a();
1043  result_type __sum = __p.a() + __p.b();
1044  while (__f != __t)
1045  {
1046  result_type __x = std::sin(this->_M_ud(__urng));
1047  *__f++ = (__x * __dif + __sum) / result_type(2);
1048  }
1049  }
1050 
1051  template<typename _RealType, typename _CharT, typename _Traits>
1054  const arcsine_distribution<_RealType>& __x)
1055  {
1056  typedef std::basic_ostream<_CharT, _Traits> __ostream_type;
1057  typedef typename __ostream_type::ios_base __ios_base;
1058 
1059  const typename __ios_base::fmtflags __flags = __os.flags();
1060  const _CharT __fill = __os.fill();
1061  const std::streamsize __precision = __os.precision();
1062  const _CharT __space = __os.widen(' ');
1064  __os.fill(__space);
1066 
1067  __os << __x.a() << __space << __x.b();
1068  __os << __space << __x._M_ud;
1069 
1070  __os.flags(__flags);
1071  __os.fill(__fill);
1072  __os.precision(__precision);
1073  return __os;
1074  }
1075 
1076  template<typename _RealType, typename _CharT, typename _Traits>
1079  arcsine_distribution<_RealType>& __x)
1080  {
1081  typedef std::basic_istream<_CharT, _Traits> __istream_type;
1082  typedef typename __istream_type::ios_base __ios_base;
1083 
1084  const typename __ios_base::fmtflags __flags = __is.flags();
1086 
1087  _RealType __a, __b;
1088  __is >> __a >> __b;
1089  __is >> __x._M_ud;
1090  __x.param(typename arcsine_distribution<_RealType>::
1091  param_type(__a, __b));
1092 
1093  __is.flags(__flags);
1094  return __is;
1095  }
1096 
1097 
1098  template<typename _RealType>
1099  template<typename _UniformRandomNumberGenerator>
1100  typename hoyt_distribution<_RealType>::result_type
1101  hoyt_distribution<_RealType>::
1102  operator()(_UniformRandomNumberGenerator& __urng)
1103  {
1104  result_type __x = this->_M_ad(__urng);
1105  result_type __y = this->_M_ed(__urng);
1106  return (result_type(2) * this->q()
1107  / (result_type(1) + this->q() * this->q()))
1108  * std::sqrt(this->omega() * __x * __y);
1109  }
1110 
1111  template<typename _RealType>
1112  template<typename _UniformRandomNumberGenerator>
1113  typename hoyt_distribution<_RealType>::result_type
1114  hoyt_distribution<_RealType>::
1115  operator()(_UniformRandomNumberGenerator& __urng,
1116  const param_type& __p)
1117  {
1118  result_type __q2 = __p.q() * __p.q();
1119  result_type __num = result_type(0.5L) * (result_type(1) + __q2);
1120  typename __gnu_cxx::arcsine_distribution<result_type>::param_type
1121  __pa(__num, __num / __q2);
1122  result_type __x = this->_M_ad(__pa, __urng);
1123  result_type __y = this->_M_ed(__urng);
1124  return (result_type(2) * __p.q() / (result_type(1) + __q2))
1125  * std::sqrt(__p.omega() * __x * __y);
1126  }
1127 
1128  template<typename _RealType>
1129  template<typename _OutputIterator,
1130  typename _UniformRandomNumberGenerator>
1131  void
1132  hoyt_distribution<_RealType>::
1133  __generate_impl(_OutputIterator __f, _OutputIterator __t,
1134  _UniformRandomNumberGenerator& __urng,
1135  const param_type& __p)
1136  {
1137  __glibcxx_function_requires(_OutputIteratorConcept<_OutputIterator,
1138  result_type>)
1139 
1140  result_type __2q = result_type(2) * __p.q();
1141  result_type __q2 = __p.q() * __p.q();
1142  result_type __q2p1 = result_type(1) + __q2;
1143  result_type __num = result_type(0.5L) * __q2p1;
1144  result_type __omega = __p.omega();
1145  typename __gnu_cxx::arcsine_distribution<result_type>::param_type
1146  __pa(__num, __num / __q2);
1147  while (__f != __t)
1148  {
1149  result_type __x = this->_M_ad(__pa, __urng);
1150  result_type __y = this->_M_ed(__urng);
1151  *__f++ = (__2q / __q2p1) * std::sqrt(__omega * __x * __y);
1152  }
1153  }
1154 
1155  template<typename _RealType, typename _CharT, typename _Traits>
1158  const hoyt_distribution<_RealType>& __x)
1159  {
1160  typedef std::basic_ostream<_CharT, _Traits> __ostream_type;
1161  typedef typename __ostream_type::ios_base __ios_base;
1162 
1163  const typename __ios_base::fmtflags __flags = __os.flags();
1164  const _CharT __fill = __os.fill();
1165  const std::streamsize __precision = __os.precision();
1166  const _CharT __space = __os.widen(' ');
1168  __os.fill(__space);
1170 
1171  __os << __x.q() << __space << __x.omega();
1172  __os << __space << __x._M_ad;
1173  __os << __space << __x._M_ed;
1174 
1175  __os.flags(__flags);
1176  __os.fill(__fill);
1177  __os.precision(__precision);
1178  return __os;
1179  }
1180 
1181  template<typename _RealType, typename _CharT, typename _Traits>
1184  hoyt_distribution<_RealType>& __x)
1185  {
1186  typedef std::basic_istream<_CharT, _Traits> __istream_type;
1187  typedef typename __istream_type::ios_base __ios_base;
1188 
1189  const typename __ios_base::fmtflags __flags = __is.flags();
1191 
1192  _RealType __q, __omega;
1193  __is >> __q >> __omega;
1194  __is >> __x._M_ad;
1195  __is >> __x._M_ed;
1196  __x.param(typename hoyt_distribution<_RealType>::
1197  param_type(__q, __omega));
1198 
1199  __is.flags(__flags);
1200  return __is;
1201  }
1202 
1203 
1204  template<typename _RealType>
1205  template<typename _OutputIterator,
1206  typename _UniformRandomNumberGenerator>
1207  void
1208  triangular_distribution<_RealType>::
1209  __generate_impl(_OutputIterator __f, _OutputIterator __t,
1210  _UniformRandomNumberGenerator& __urng,
1211  const param_type& __param)
1212  {
1213  __glibcxx_function_requires(_OutputIteratorConcept<_OutputIterator,
1214  result_type>)
1215 
1216  while (__f != __t)
1217  *__f++ = this->operator()(__urng, __param);
1218  }
1219 
1220  template<typename _RealType, typename _CharT, typename _Traits>
1221  std::basic_ostream<_CharT, _Traits>&
1222  operator<<(std::basic_ostream<_CharT, _Traits>& __os,
1223  const __gnu_cxx::triangular_distribution<_RealType>& __x)
1224  {
1225  typedef std::basic_ostream<_CharT, _Traits> __ostream_type;
1226  typedef typename __ostream_type::ios_base __ios_base;
1227 
1228  const typename __ios_base::fmtflags __flags = __os.flags();
1229  const _CharT __fill = __os.fill();
1230  const std::streamsize __precision = __os.precision();
1231  const _CharT __space = __os.widen(' ');
1233  __os.fill(__space);
1235 
1236  __os << __x.a() << __space << __x.b() << __space << __x.c();
1237 
1238  __os.flags(__flags);
1239  __os.fill(__fill);
1240  __os.precision(__precision);
1241  return __os;
1242  }
1243 
1244  template<typename _RealType, typename _CharT, typename _Traits>
1247  __gnu_cxx::triangular_distribution<_RealType>& __x)
1248  {
1249  typedef std::basic_istream<_CharT, _Traits> __istream_type;
1250  typedef typename __istream_type::ios_base __ios_base;
1251 
1252  const typename __ios_base::fmtflags __flags = __is.flags();
1254 
1255  _RealType __a, __b, __c;
1256  __is >> __a >> __b >> __c;
1257  __x.param(typename __gnu_cxx::triangular_distribution<_RealType>::
1258  param_type(__a, __b, __c));
1259 
1260  __is.flags(__flags);
1261  return __is;
1262  }
1263 
1264 
1265  template<typename _RealType>
1266  template<typename _UniformRandomNumberGenerator>
1267  typename von_mises_distribution<_RealType>::result_type
1268  von_mises_distribution<_RealType>::
1269  operator()(_UniformRandomNumberGenerator& __urng,
1270  const param_type& __p)
1271  {
1272  const result_type __pi
1273  = __gnu_cxx::__math_constants<result_type>::__pi;
1274  std::__detail::_Adaptor<_UniformRandomNumberGenerator, result_type>
1275  __aurng(__urng);
1276 
1277  result_type __f;
1278  while (1)
1279  {
1280  result_type __rnd = std::cos(__pi * __aurng());
1281  __f = (result_type(1) + __p._M_r * __rnd) / (__p._M_r + __rnd);
1282  result_type __c = __p._M_kappa * (__p._M_r - __f);
1283 
1284  result_type __rnd2 = __aurng();
1285  if (__c * (result_type(2) - __c) > __rnd2)
1286  break;
1287  if (std::log(__c / __rnd2) >= __c - result_type(1))
1288  break;
1289  }
1290 
1291  result_type __res = std::acos(__f);
1292 #if _GLIBCXX_USE_C99_MATH_FUNCS
1293  __res = std::copysign(__res, __aurng() - result_type(0.5));
1294 #else
1295  if (__aurng() < result_type(0.5))
1296  __res = -__res;
1297 #endif
1298  __res += __p._M_mu;
1299  if (__res > __pi)
1300  __res -= result_type(2) * __pi;
1301  else if (__res < -__pi)
1302  __res += result_type(2) * __pi;
1303  return __res;
1304  }
1305 
1306  template<typename _RealType>
1307  template<typename _OutputIterator,
1308  typename _UniformRandomNumberGenerator>
1309  void
1310  von_mises_distribution<_RealType>::
1311  __generate_impl(_OutputIterator __f, _OutputIterator __t,
1312  _UniformRandomNumberGenerator& __urng,
1313  const param_type& __param)
1314  {
1315  __glibcxx_function_requires(_OutputIteratorConcept<_OutputIterator,
1316  result_type>)
1317 
1318  while (__f != __t)
1319  *__f++ = this->operator()(__urng, __param);
1320  }
1321 
1322  template<typename _RealType, typename _CharT, typename _Traits>
1323  std::basic_ostream<_CharT, _Traits>&
1324  operator<<(std::basic_ostream<_CharT, _Traits>& __os,
1325  const __gnu_cxx::von_mises_distribution<_RealType>& __x)
1326  {
1327  typedef std::basic_ostream<_CharT, _Traits> __ostream_type;
1328  typedef typename __ostream_type::ios_base __ios_base;
1329 
1330  const typename __ios_base::fmtflags __flags = __os.flags();
1331  const _CharT __fill = __os.fill();
1332  const std::streamsize __precision = __os.precision();
1333  const _CharT __space = __os.widen(' ');
1335  __os.fill(__space);
1337 
1338  __os << __x.mu() << __space << __x.kappa();
1339 
1340  __os.flags(__flags);
1341  __os.fill(__fill);
1342  __os.precision(__precision);
1343  return __os;
1344  }
1345 
1346  template<typename _RealType, typename _CharT, typename _Traits>
1349  __gnu_cxx::von_mises_distribution<_RealType>& __x)
1350  {
1351  typedef std::basic_istream<_CharT, _Traits> __istream_type;
1352  typedef typename __istream_type::ios_base __ios_base;
1353 
1354  const typename __ios_base::fmtflags __flags = __is.flags();
1356 
1357  _RealType __mu, __kappa;
1358  __is >> __mu >> __kappa;
1359  __x.param(typename __gnu_cxx::von_mises_distribution<_RealType>::
1360  param_type(__mu, __kappa));
1361 
1362  __is.flags(__flags);
1363  return __is;
1364  }
1365 
1366 
1367  template<typename _UIntType>
1368  template<typename _UniformRandomNumberGenerator>
1369  typename hypergeometric_distribution<_UIntType>::result_type
1370  hypergeometric_distribution<_UIntType>::
1371  operator()(_UniformRandomNumberGenerator& __urng,
1372  const param_type& __param)
1373  {
1374  std::__detail::_Adaptor<_UniformRandomNumberGenerator, double>
1375  __aurng(__urng);
1376 
1377  result_type __a = __param.successful_size();
1378  result_type __b = __param.total_size();
1379  result_type __k = 0;
1380 
1381  if (__param.total_draws() < __param.total_size() / 2)
1382  {
1383  for (result_type __i = 0; __i < __param.total_draws(); ++__i)
1384  {
1385  if (__b * __aurng() < __a)
1386  {
1387  ++__k;
1388  if (__k == __param.successful_size())
1389  return __k;
1390  --__a;
1391  }
1392  --__b;
1393  }
1394  return __k;
1395  }
1396  else
1397  {
1398  for (result_type __i = 0; __i < __param.unsuccessful_size(); ++__i)
1399  {
1400  if (__b * __aurng() < __a)
1401  {
1402  ++__k;
1403  if (__k == __param.successful_size())
1404  return __param.successful_size() - __k;
1405  --__a;
1406  }
1407  --__b;
1408  }
1409  return __param.successful_size() - __k;
1410  }
1411  }
1412 
1413  template<typename _UIntType>
1414  template<typename _OutputIterator,
1415  typename _UniformRandomNumberGenerator>
1416  void
1417  hypergeometric_distribution<_UIntType>::
1418  __generate_impl(_OutputIterator __f, _OutputIterator __t,
1419  _UniformRandomNumberGenerator& __urng,
1420  const param_type& __param)
1421  {
1422  __glibcxx_function_requires(_OutputIteratorConcept<_OutputIterator,
1423  result_type>)
1424 
1425  while (__f != __t)
1426  *__f++ = this->operator()(__urng);
1427  }
1428 
1429  template<typename _UIntType, typename _CharT, typename _Traits>
1430  std::basic_ostream<_CharT, _Traits>&
1431  operator<<(std::basic_ostream<_CharT, _Traits>& __os,
1432  const __gnu_cxx::hypergeometric_distribution<_UIntType>& __x)
1433  {
1434  typedef std::basic_ostream<_CharT, _Traits> __ostream_type;
1435  typedef typename __ostream_type::ios_base __ios_base;
1436 
1437  const typename __ios_base::fmtflags __flags = __os.flags();
1438  const _CharT __fill = __os.fill();
1439  const std::streamsize __precision = __os.precision();
1440  const _CharT __space = __os.widen(' ');
1442  __os.fill(__space);
1444 
1445  __os << __x.total_size() << __space << __x.successful_size() << __space
1446  << __x.total_draws();
1447 
1448  __os.flags(__flags);
1449  __os.fill(__fill);
1450  __os.precision(__precision);
1451  return __os;
1452  }
1453 
1454  template<typename _UIntType, typename _CharT, typename _Traits>
1457  __gnu_cxx::hypergeometric_distribution<_UIntType>& __x)
1458  {
1459  typedef std::basic_istream<_CharT, _Traits> __istream_type;
1460  typedef typename __istream_type::ios_base __ios_base;
1461 
1462  const typename __ios_base::fmtflags __flags = __is.flags();
1464 
1465  _UIntType __total_size, __successful_size, __total_draws;
1466  __is >> __total_size >> __successful_size >> __total_draws;
1467  __x.param(typename __gnu_cxx::hypergeometric_distribution<_UIntType>::
1468  param_type(__total_size, __successful_size, __total_draws));
1469 
1470  __is.flags(__flags);
1471  return __is;
1472  }
1473 
1474 
1475  template<typename _RealType>
1476  template<typename _UniformRandomNumberGenerator>
1477  typename logistic_distribution<_RealType>::result_type
1478  logistic_distribution<_RealType>::
1479  operator()(_UniformRandomNumberGenerator& __urng,
1480  const param_type& __p)
1481  {
1482  std::__detail::_Adaptor<_UniformRandomNumberGenerator, result_type>
1483  __aurng(__urng);
1484 
1485  result_type __arg = result_type(1);
1486  while (__arg == result_type(1) || __arg == result_type(0))
1487  __arg = __aurng();
1488  return __p.a()
1489  + __p.b() * std::log(__arg / (result_type(1) - __arg));
1490  }
1491 
1492  template<typename _RealType>
1493  template<typename _OutputIterator,
1494  typename _UniformRandomNumberGenerator>
1495  void
1496  logistic_distribution<_RealType>::
1497  __generate_impl(_OutputIterator __f, _OutputIterator __t,
1498  _UniformRandomNumberGenerator& __urng,
1499  const param_type& __p)
1500  {
1501  __glibcxx_function_requires(_OutputIteratorConcept<_OutputIterator,
1502  result_type>)
1503 
1504  std::__detail::_Adaptor<_UniformRandomNumberGenerator, result_type>
1505  __aurng(__urng);
1506 
1507  while (__f != __t)
1508  {
1509  result_type __arg = result_type(1);
1510  while (__arg == result_type(1) || __arg == result_type(0))
1511  __arg = __aurng();
1512  *__f++ = __p.a()
1513  + __p.b() * std::log(__arg / (result_type(1) - __arg));
1514  }
1515  }
1516 
1517  template<typename _RealType, typename _CharT, typename _Traits>
1520  const logistic_distribution<_RealType>& __x)
1521  {
1522  typedef std::basic_ostream<_CharT, _Traits> __ostream_type;
1523  typedef typename __ostream_type::ios_base __ios_base;
1524 
1525  const typename __ios_base::fmtflags __flags = __os.flags();
1526  const _CharT __fill = __os.fill();
1527  const std::streamsize __precision = __os.precision();
1528  const _CharT __space = __os.widen(' ');
1530  __os.fill(__space);
1532 
1533  __os << __x.a() << __space << __x.b();
1534 
1535  __os.flags(__flags);
1536  __os.fill(__fill);
1537  __os.precision(__precision);
1538  return __os;
1539  }
1540 
1541  template<typename _RealType, typename _CharT, typename _Traits>
1544  logistic_distribution<_RealType>& __x)
1545  {
1546  typedef std::basic_istream<_CharT, _Traits> __istream_type;
1547  typedef typename __istream_type::ios_base __ios_base;
1548 
1549  const typename __ios_base::fmtflags __flags = __is.flags();
1551 
1552  _RealType __a, __b;
1553  __is >> __a >> __b;
1554  __x.param(typename logistic_distribution<_RealType>::
1555  param_type(__a, __b));
1556 
1557  __is.flags(__flags);
1558  return __is;
1559  }
1560 
1561 
1562  namespace {
1563 
1564  // Helper class for the uniform_on_sphere_distribution generation
1565  // function.
1566  template<std::size_t _Dimen, typename _RealType>
1567  class uniform_on_sphere_helper
1568  {
1569  typedef typename uniform_on_sphere_distribution<_Dimen, _RealType>::
1570  result_type result_type;
1571 
1572  public:
1573  template<typename _NormalDistribution,
1574  typename _UniformRandomNumberGenerator>
1575  result_type operator()(_NormalDistribution& __nd,
1576  _UniformRandomNumberGenerator& __urng)
1577  {
1578  result_type __ret;
1579  typename result_type::value_type __norm;
1580 
1581  do
1582  {
1583  auto __sum = _RealType(0);
1584 
1585  std::generate(__ret.begin(), __ret.end(),
1586  [&__nd, &__urng, &__sum](){
1587  _RealType __t = __nd(__urng);
1588  __sum += __t * __t;
1589  return __t; });
1590  __norm = std::sqrt(__sum);
1591  }
1592  while (__norm == _RealType(0) || ! __builtin_isfinite(__norm));
1593 
1594  std::transform(__ret.begin(), __ret.end(), __ret.begin(),
1595  [__norm](_RealType __val){ return __val / __norm; });
1596 
1597  return __ret;
1598  }
1599  };
1600 
1601 
1602  template<typename _RealType>
1603  class uniform_on_sphere_helper<2, _RealType>
1604  {
1605  typedef typename uniform_on_sphere_distribution<2, _RealType>::
1606  result_type result_type;
1607 
1608  public:
1609  template<typename _NormalDistribution,
1610  typename _UniformRandomNumberGenerator>
1611  result_type operator()(_NormalDistribution&,
1612  _UniformRandomNumberGenerator& __urng)
1613  {
1614  result_type __ret;
1615  _RealType __sq;
1616  std::__detail::_Adaptor<_UniformRandomNumberGenerator,
1617  _RealType> __aurng(__urng);
1618 
1619  do
1620  {
1621  __ret[0] = _RealType(2) * __aurng() - _RealType(1);
1622  __ret[1] = _RealType(2) * __aurng() - _RealType(1);
1623 
1624  __sq = __ret[0] * __ret[0] + __ret[1] * __ret[1];
1625  }
1626  while (__sq == _RealType(0) || __sq > _RealType(1));
1627 
1628 #if _GLIBCXX_USE_C99_MATH_FUNCS
1629  // Yes, we do not just use sqrt(__sq) because hypot() is more
1630  // accurate.
1631  auto __norm = std::hypot(__ret[0], __ret[1]);
1632 #else
1633  auto __norm = std::sqrt(__sq);
1634 #endif
1635  __ret[0] /= __norm;
1636  __ret[1] /= __norm;
1637 
1638  return __ret;
1639  }
1640  };
1641 
1642  }
1643 
1644 
1645  template<std::size_t _Dimen, typename _RealType>
1646  template<typename _UniformRandomNumberGenerator>
1647  typename uniform_on_sphere_distribution<_Dimen, _RealType>::result_type
1648  uniform_on_sphere_distribution<_Dimen, _RealType>::
1649  operator()(_UniformRandomNumberGenerator& __urng,
1650  const param_type& __p)
1651  {
1652  uniform_on_sphere_helper<_Dimen, _RealType> __helper;
1653  return __helper(_M_nd, __urng);
1654  }
1655 
1656  template<std::size_t _Dimen, typename _RealType>
1657  template<typename _OutputIterator,
1658  typename _UniformRandomNumberGenerator>
1659  void
1660  uniform_on_sphere_distribution<_Dimen, _RealType>::
1661  __generate_impl(_OutputIterator __f, _OutputIterator __t,
1662  _UniformRandomNumberGenerator& __urng,
1663  const param_type& __param)
1664  {
1665  __glibcxx_function_requires(_OutputIteratorConcept<_OutputIterator,
1666  result_type>)
1667 
1668  while (__f != __t)
1669  *__f++ = this->operator()(__urng, __param);
1670  }
1671 
1672  template<std::size_t _Dimen, typename _RealType, typename _CharT,
1673  typename _Traits>
1674  std::basic_ostream<_CharT, _Traits>&
1675  operator<<(std::basic_ostream<_CharT, _Traits>& __os,
1676  const __gnu_cxx::uniform_on_sphere_distribution<_Dimen,
1677  _RealType>& __x)
1678  {
1679  return __os << __x._M_nd;
1680  }
1681 
1682  template<std::size_t _Dimen, typename _RealType, typename _CharT,
1683  typename _Traits>
1686  __gnu_cxx::uniform_on_sphere_distribution<_Dimen,
1687  _RealType>& __x)
1688  {
1689  return __is >> __x._M_nd;
1690  }
1691 
1692 
1693  namespace {
1694 
1695  // Helper class for the uniform_inside_sphere_distribution generation
1696  // function.
1697  template<std::size_t _Dimen, bool _SmallDimen, typename _RealType>
1698  class uniform_inside_sphere_helper;
1699 
1700  template<std::size_t _Dimen, typename _RealType>
1701  class uniform_inside_sphere_helper<_Dimen, false, _RealType>
1702  {
1703  using result_type
1704  = typename uniform_inside_sphere_distribution<_Dimen, _RealType>::
1705  result_type;
1706 
1707  public:
1708  template<typename _UniformOnSphereDistribution,
1709  typename _UniformRandomNumberGenerator>
1710  result_type
1711  operator()(_UniformOnSphereDistribution& __uosd,
1712  _UniformRandomNumberGenerator& __urng,
1713  _RealType __radius)
1714  {
1715  std::__detail::_Adaptor<_UniformRandomNumberGenerator,
1716  _RealType> __aurng(__urng);
1717 
1718  _RealType __pow = 1 / _RealType(_Dimen);
1719  _RealType __urt = __radius * std::pow(__aurng(), __pow);
1720  result_type __ret = __uosd(__aurng);
1721 
1722  std::transform(__ret.begin(), __ret.end(), __ret.begin(),
1723  [__urt](_RealType __val)
1724  { return __val * __urt; });
1725 
1726  return __ret;
1727  }
1728  };
1729 
1730  // Helper class for the uniform_inside_sphere_distribution generation
1731  // function specialized for small dimensions.
1732  template<std::size_t _Dimen, typename _RealType>
1733  class uniform_inside_sphere_helper<_Dimen, true, _RealType>
1734  {
1735  using result_type
1736  = typename uniform_inside_sphere_distribution<_Dimen, _RealType>::
1737  result_type;
1738 
1739  public:
1740  template<typename _UniformOnSphereDistribution,
1741  typename _UniformRandomNumberGenerator>
1742  result_type
1743  operator()(_UniformOnSphereDistribution&,
1744  _UniformRandomNumberGenerator& __urng,
1745  _RealType __radius)
1746  {
1747  result_type __ret;
1748  _RealType __sq;
1749  _RealType __radsq = __radius * __radius;
1750  std::__detail::_Adaptor<_UniformRandomNumberGenerator,
1751  _RealType> __aurng(__urng);
1752 
1753  do
1754  {
1755  __sq = _RealType(0);
1756  for (int i = 0; i < _Dimen; ++i)
1757  {
1758  __ret[i] = _RealType(2) * __aurng() - _RealType(1);
1759  __sq += __ret[i] * __ret[i];
1760  }
1761  }
1762  while (__sq > _RealType(1));
1763 
1764  for (int i = 0; i < _Dimen; ++i)
1765  __ret[i] *= __radius;
1766 
1767  return __ret;
1768  }
1769  };
1770  } // namespace
1771 
1772  //
1773  // Experiments have shown that rejection is more efficient than transform
1774  // for dimensions less than 8.
1775  //
1776  template<std::size_t _Dimen, typename _RealType>
1777  template<typename _UniformRandomNumberGenerator>
1778  typename uniform_inside_sphere_distribution<_Dimen, _RealType>::result_type
1779  uniform_inside_sphere_distribution<_Dimen, _RealType>::
1780  operator()(_UniformRandomNumberGenerator& __urng,
1781  const param_type& __p)
1782  {
1783  uniform_inside_sphere_helper<_Dimen, _Dimen < 8, _RealType> __helper;
1784  return __helper(_M_uosd, __urng, __p.radius());
1785  }
1786 
1787  template<std::size_t _Dimen, typename _RealType>
1788  template<typename _OutputIterator,
1789  typename _UniformRandomNumberGenerator>
1790  void
1791  uniform_inside_sphere_distribution<_Dimen, _RealType>::
1792  __generate_impl(_OutputIterator __f, _OutputIterator __t,
1793  _UniformRandomNumberGenerator& __urng,
1794  const param_type& __param)
1795  {
1796  __glibcxx_function_requires(_OutputIteratorConcept<_OutputIterator,
1797  result_type>)
1798 
1799  while (__f != __t)
1800  *__f++ = this->operator()(__urng, __param);
1801  }
1802 
1803  template<std::size_t _Dimen, typename _RealType, typename _CharT,
1804  typename _Traits>
1805  std::basic_ostream<_CharT, _Traits>&
1806  operator<<(std::basic_ostream<_CharT, _Traits>& __os,
1807  const __gnu_cxx::uniform_inside_sphere_distribution<_Dimen,
1808  _RealType>& __x)
1809  {
1810  typedef std::basic_ostream<_CharT, _Traits> __ostream_type;
1811  typedef typename __ostream_type::ios_base __ios_base;
1812 
1813  const typename __ios_base::fmtflags __flags = __os.flags();
1814  const _CharT __fill = __os.fill();
1815  const std::streamsize __precision = __os.precision();
1816  const _CharT __space = __os.widen(' ');
1818  __os.fill(__space);
1820 
1821  __os << __x.radius() << __space << __x._M_uosd;
1822 
1823  __os.flags(__flags);
1824  __os.fill(__fill);
1825  __os.precision(__precision);
1826 
1827  return __os;
1828  }
1829 
1830  template<std::size_t _Dimen, typename _RealType, typename _CharT,
1831  typename _Traits>
1834  __gnu_cxx::uniform_inside_sphere_distribution<_Dimen,
1835  _RealType>& __x)
1836  {
1837  typedef std::basic_istream<_CharT, _Traits> __istream_type;
1838  typedef typename __istream_type::ios_base __ios_base;
1839 
1840  const typename __ios_base::fmtflags __flags = __is.flags();
1842 
1843  _RealType __radius_val;
1844  __is >> __radius_val >> __x._M_uosd;
1845  __x.param(typename uniform_inside_sphere_distribution<_Dimen, _RealType>::
1846  param_type(__radius_val));
1847 
1848  __is.flags(__flags);
1849 
1850  return __is;
1851  }
1852 
1853 _GLIBCXX_END_NAMESPACE_VERSION
1854 } // namespace __gnu_cxx
1855 
1856 
1857 #endif // _EXT_RANDOM_TCC
complex< _Tp > sin(const complex< _Tp > &)
Return complex sine of z.
Definition: complex:1197
complex< _Tp > log(const complex< _Tp > &)
Return complex natural logarithm of z.
Definition: complex:1162
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 > cos(const complex< _Tp > &)
Return complex cosine of z.
Definition: complex:1079
complex< _Tp > sqrt(const complex< _Tp > &)
Return complex square root of z.
Definition: complex:1271
ISO C++ entities toplevel namespace is std.
std::basic_istream< _CharT, _Traits > & operator>>(std::basic_istream< _CharT, _Traits > &__is, bitset< _Nb > &__x)
Global I/O operators for bitsets.
Definition: bitset:1602
ptrdiff_t streamsize
Integral type for I/O operation counts and buffer sizes.
Definition: postypes.h:73
std::basic_ostream< _CharT, _Traits > & operator<<(std::basic_ostream< _CharT, _Traits > &__os, const bitset< _Nb > &__x)
Global I/O operators for bitsets.
Definition: bitset:1692
ios_base & scientific(ios_base &__base)
Calls base.setf(ios_base::scientific, ios_base::floatfield).
Definition: ios_base.h:1127
ios_base & dec(ios_base &__base)
Calls base.setf(ios_base::dec, ios_base::basefield).
Definition: ios_base.h:1094
ios_base & left(ios_base &__base)
Calls base.setf(ios_base::left, ios_base::adjustfield).
Definition: ios_base.h:1077
constexpr void advance(_InputIterator &__i, _Distance __n)
A generalization of pointer arithmetic.
ios_base & skipws(ios_base &__base)
Calls base.setf(ios_base::skipws).
Definition: ios_base.h:1020
ios_base & fixed(ios_base &__base)
Calls base.setf(ios_base::fixed, ios_base::floatfield).
Definition: ios_base.h:1119
GNU extensions for public use.
A standard container for storing a fixed size sequence of elements.
Definition: array:103
Template class basic_istream.
Definition: istream:63
Template class basic_ostream.
Definition: ostream.h:67
Properties of fundamental types.
Definition: limits:320
fmtflags flags() const
Access to format flags.
Definition: ios_base.h:694