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