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