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