libstdc++
ext/random.tcc
Go to the documentation of this file.
1 // Random number extensions -*- C++ -*-
2 
3 // Copyright (C) 2012-2021 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>
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>
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_lower"));
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>
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  // The param_type temporary is built with a private constructor,
713  // to skip the Cholesky decomposition that would be performed
714  // otherwise.
715  __x.param(typename normal_mv_distribution<_Dimen, _RealType>::
716  param_type(__mean, __varcov));
717 
718  __is.flags(__flags);
719  return __is;
720  }
721 
722 
723  template<typename _RealType>
724  template<typename _OutputIterator,
725  typename _UniformRandomNumberGenerator>
726  void
727  rice_distribution<_RealType>::
728  __generate_impl(_OutputIterator __f, _OutputIterator __t,
729  _UniformRandomNumberGenerator& __urng,
730  const param_type& __p)
731  {
732  __glibcxx_function_requires(_OutputIteratorConcept<_OutputIterator,
733  result_type>)
734 
735  while (__f != __t)
736  {
738  __px(__p.nu(), __p.sigma()), __py(result_type(0), __p.sigma());
739  result_type __x = this->_M_ndx(__px, __urng);
740  result_type __y = this->_M_ndy(__py, __urng);
741 #if _GLIBCXX_USE_C99_MATH_TR1
742  *__f++ = std::hypot(__x, __y);
743 #else
744  *__f++ = std::sqrt(__x * __x + __y * __y);
745 #endif
746  }
747  }
748 
749  template<typename _RealType, typename _CharT, typename _Traits>
752  const rice_distribution<_RealType>& __x)
753  {
754  typedef std::basic_ostream<_CharT, _Traits> __ostream_type;
755  typedef typename __ostream_type::ios_base __ios_base;
756 
757  const typename __ios_base::fmtflags __flags = __os.flags();
758  const _CharT __fill = __os.fill();
759  const std::streamsize __precision = __os.precision();
760  const _CharT __space = __os.widen(' ');
762  __os.fill(__space);
764 
765  __os << __x.nu() << __space << __x.sigma();
766  __os << __space << __x._M_ndx;
767  __os << __space << __x._M_ndy;
768 
769  __os.flags(__flags);
770  __os.fill(__fill);
771  __os.precision(__precision);
772  return __os;
773  }
774 
775  template<typename _RealType, typename _CharT, typename _Traits>
778  rice_distribution<_RealType>& __x)
779  {
780  typedef std::basic_istream<_CharT, _Traits> __istream_type;
781  typedef typename __istream_type::ios_base __ios_base;
782 
783  const typename __ios_base::fmtflags __flags = __is.flags();
785 
786  _RealType __nu_val, __sigma_val;
787  __is >> __nu_val >> __sigma_val;
788  __is >> __x._M_ndx;
789  __is >> __x._M_ndy;
790  __x.param(typename rice_distribution<_RealType>::
791  param_type(__nu_val, __sigma_val));
792 
793  __is.flags(__flags);
794  return __is;
795  }
796 
797 
798  template<typename _RealType>
799  template<typename _OutputIterator,
800  typename _UniformRandomNumberGenerator>
801  void
802  nakagami_distribution<_RealType>::
803  __generate_impl(_OutputIterator __f, _OutputIterator __t,
804  _UniformRandomNumberGenerator& __urng,
805  const param_type& __p)
806  {
807  __glibcxx_function_requires(_OutputIteratorConcept<_OutputIterator,
808  result_type>)
809 
810  typename std::gamma_distribution<result_type>::param_type
811  __pg(__p.mu(), __p.omega() / __p.mu());
812  while (__f != __t)
813  *__f++ = std::sqrt(this->_M_gd(__pg, __urng));
814  }
815 
816  template<typename _RealType, typename _CharT, typename _Traits>
817  std::basic_ostream<_CharT, _Traits>&
818  operator<<(std::basic_ostream<_CharT, _Traits>& __os,
819  const nakagami_distribution<_RealType>& __x)
820  {
821  typedef std::basic_ostream<_CharT, _Traits> __ostream_type;
822  typedef typename __ostream_type::ios_base __ios_base;
823 
824  const typename __ios_base::fmtflags __flags = __os.flags();
825  const _CharT __fill = __os.fill();
826  const std::streamsize __precision = __os.precision();
827  const _CharT __space = __os.widen(' ');
829  __os.fill(__space);
831 
832  __os << __x.mu() << __space << __x.omega();
833  __os << __space << __x._M_gd;
834 
835  __os.flags(__flags);
836  __os.fill(__fill);
837  __os.precision(__precision);
838  return __os;
839  }
840 
841  template<typename _RealType, typename _CharT, typename _Traits>
844  nakagami_distribution<_RealType>& __x)
845  {
846  typedef std::basic_istream<_CharT, _Traits> __istream_type;
847  typedef typename __istream_type::ios_base __ios_base;
848 
849  const typename __ios_base::fmtflags __flags = __is.flags();
851 
852  _RealType __mu_val, __omega_val;
853  __is >> __mu_val >> __omega_val;
854  __is >> __x._M_gd;
855  __x.param(typename nakagami_distribution<_RealType>::
856  param_type(__mu_val, __omega_val));
857 
858  __is.flags(__flags);
859  return __is;
860  }
861 
862 
863  template<typename _RealType>
864  template<typename _OutputIterator,
865  typename _UniformRandomNumberGenerator>
866  void
867  pareto_distribution<_RealType>::
868  __generate_impl(_OutputIterator __f, _OutputIterator __t,
869  _UniformRandomNumberGenerator& __urng,
870  const param_type& __p)
871  {
872  __glibcxx_function_requires(_OutputIteratorConcept<_OutputIterator,
873  result_type>)
874 
875  result_type __mu_val = __p.mu();
876  result_type __malphinv = -result_type(1) / __p.alpha();
877  while (__f != __t)
878  *__f++ = __mu_val * std::pow(this->_M_ud(__urng), __malphinv);
879  }
880 
881  template<typename _RealType, typename _CharT, typename _Traits>
882  std::basic_ostream<_CharT, _Traits>&
883  operator<<(std::basic_ostream<_CharT, _Traits>& __os,
884  const pareto_distribution<_RealType>& __x)
885  {
886  typedef std::basic_ostream<_CharT, _Traits> __ostream_type;
887  typedef typename __ostream_type::ios_base __ios_base;
888 
889  const typename __ios_base::fmtflags __flags = __os.flags();
890  const _CharT __fill = __os.fill();
891  const std::streamsize __precision = __os.precision();
892  const _CharT __space = __os.widen(' ');
894  __os.fill(__space);
896 
897  __os << __x.alpha() << __space << __x.mu();
898  __os << __space << __x._M_ud;
899 
900  __os.flags(__flags);
901  __os.fill(__fill);
902  __os.precision(__precision);
903  return __os;
904  }
905 
906  template<typename _RealType, typename _CharT, typename _Traits>
909  pareto_distribution<_RealType>& __x)
910  {
911  typedef std::basic_istream<_CharT, _Traits> __istream_type;
912  typedef typename __istream_type::ios_base __ios_base;
913 
914  const typename __ios_base::fmtflags __flags = __is.flags();
916 
917  _RealType __alpha_val, __mu_val;
918  __is >> __alpha_val >> __mu_val;
919  __is >> __x._M_ud;
920  __x.param(typename pareto_distribution<_RealType>::
921  param_type(__alpha_val, __mu_val));
922 
923  __is.flags(__flags);
924  return __is;
925  }
926 
927 
928  template<typename _RealType>
929  template<typename _UniformRandomNumberGenerator>
930  typename k_distribution<_RealType>::result_type
931  k_distribution<_RealType>::
932  operator()(_UniformRandomNumberGenerator& __urng)
933  {
934  result_type __x = this->_M_gd1(__urng);
935  result_type __y = this->_M_gd2(__urng);
936  return std::sqrt(__x * __y);
937  }
938 
939  template<typename _RealType>
940  template<typename _UniformRandomNumberGenerator>
941  typename k_distribution<_RealType>::result_type
942  k_distribution<_RealType>::
943  operator()(_UniformRandomNumberGenerator& __urng,
944  const param_type& __p)
945  {
947  __p1(__p.lambda(), result_type(1) / __p.lambda()),
948  __p2(__p.nu(), __p.mu() / __p.nu());
949  result_type __x = this->_M_gd1(__p1, __urng);
950  result_type __y = this->_M_gd2(__p2, __urng);
951  return std::sqrt(__x * __y);
952  }
953 
954  template<typename _RealType>
955  template<typename _OutputIterator,
956  typename _UniformRandomNumberGenerator>
957  void
958  k_distribution<_RealType>::
959  __generate_impl(_OutputIterator __f, _OutputIterator __t,
960  _UniformRandomNumberGenerator& __urng,
961  const param_type& __p)
962  {
963  __glibcxx_function_requires(_OutputIteratorConcept<_OutputIterator,
964  result_type>)
965 
966  typename std::gamma_distribution<result_type>::param_type
967  __p1(__p.lambda(), result_type(1) / __p.lambda()),
968  __p2(__p.nu(), __p.mu() / __p.nu());
969  while (__f != __t)
970  {
971  result_type __x = this->_M_gd1(__p1, __urng);
972  result_type __y = this->_M_gd2(__p2, __urng);
973  *__f++ = std::sqrt(__x * __y);
974  }
975  }
976 
977  template<typename _RealType, typename _CharT, typename _Traits>
980  const k_distribution<_RealType>& __x)
981  {
982  typedef std::basic_ostream<_CharT, _Traits> __ostream_type;
983  typedef typename __ostream_type::ios_base __ios_base;
984 
985  const typename __ios_base::fmtflags __flags = __os.flags();
986  const _CharT __fill = __os.fill();
987  const std::streamsize __precision = __os.precision();
988  const _CharT __space = __os.widen(' ');
990  __os.fill(__space);
992 
993  __os << __x.lambda() << __space << __x.mu() << __space << __x.nu();
994  __os << __space << __x._M_gd1;
995  __os << __space << __x._M_gd2;
996 
997  __os.flags(__flags);
998  __os.fill(__fill);
999  __os.precision(__precision);
1000  return __os;
1001  }
1002 
1003  template<typename _RealType, typename _CharT, typename _Traits>
1006  k_distribution<_RealType>& __x)
1007  {
1008  typedef std::basic_istream<_CharT, _Traits> __istream_type;
1009  typedef typename __istream_type::ios_base __ios_base;
1010 
1011  const typename __ios_base::fmtflags __flags = __is.flags();
1013 
1014  _RealType __lambda_val, __mu_val, __nu_val;
1015  __is >> __lambda_val >> __mu_val >> __nu_val;
1016  __is >> __x._M_gd1;
1017  __is >> __x._M_gd2;
1018  __x.param(typename k_distribution<_RealType>::
1019  param_type(__lambda_val, __mu_val, __nu_val));
1020 
1021  __is.flags(__flags);
1022  return __is;
1023  }
1024 
1025 
1026  template<typename _RealType>
1027  template<typename _OutputIterator,
1028  typename _UniformRandomNumberGenerator>
1029  void
1030  arcsine_distribution<_RealType>::
1031  __generate_impl(_OutputIterator __f, _OutputIterator __t,
1032  _UniformRandomNumberGenerator& __urng,
1033  const param_type& __p)
1034  {
1035  __glibcxx_function_requires(_OutputIteratorConcept<_OutputIterator,
1036  result_type>)
1037 
1038  result_type __dif = __p.b() - __p.a();
1039  result_type __sum = __p.a() + __p.b();
1040  while (__f != __t)
1041  {
1042  result_type __x = std::sin(this->_M_ud(__urng));
1043  *__f++ = (__x * __dif + __sum) / result_type(2);
1044  }
1045  }
1046 
1047  template<typename _RealType, typename _CharT, typename _Traits>
1050  const arcsine_distribution<_RealType>& __x)
1051  {
1052  typedef std::basic_ostream<_CharT, _Traits> __ostream_type;
1053  typedef typename __ostream_type::ios_base __ios_base;
1054 
1055  const typename __ios_base::fmtflags __flags = __os.flags();
1056  const _CharT __fill = __os.fill();
1057  const std::streamsize __precision = __os.precision();
1058  const _CharT __space = __os.widen(' ');
1060  __os.fill(__space);
1062 
1063  __os << __x.a() << __space << __x.b();
1064  __os << __space << __x._M_ud;
1065 
1066  __os.flags(__flags);
1067  __os.fill(__fill);
1068  __os.precision(__precision);
1069  return __os;
1070  }
1071 
1072  template<typename _RealType, typename _CharT, typename _Traits>
1075  arcsine_distribution<_RealType>& __x)
1076  {
1077  typedef std::basic_istream<_CharT, _Traits> __istream_type;
1078  typedef typename __istream_type::ios_base __ios_base;
1079 
1080  const typename __ios_base::fmtflags __flags = __is.flags();
1082 
1083  _RealType __a, __b;
1084  __is >> __a >> __b;
1085  __is >> __x._M_ud;
1086  __x.param(typename arcsine_distribution<_RealType>::
1087  param_type(__a, __b));
1088 
1089  __is.flags(__flags);
1090  return __is;
1091  }
1092 
1093 
1094  template<typename _RealType>
1095  template<typename _UniformRandomNumberGenerator>
1096  typename hoyt_distribution<_RealType>::result_type
1097  hoyt_distribution<_RealType>::
1098  operator()(_UniformRandomNumberGenerator& __urng)
1099  {
1100  result_type __x = this->_M_ad(__urng);
1101  result_type __y = this->_M_ed(__urng);
1102  return (result_type(2) * this->q()
1103  / (result_type(1) + this->q() * this->q()))
1104  * std::sqrt(this->omega() * __x * __y);
1105  }
1106 
1107  template<typename _RealType>
1108  template<typename _UniformRandomNumberGenerator>
1109  typename hoyt_distribution<_RealType>::result_type
1110  hoyt_distribution<_RealType>::
1111  operator()(_UniformRandomNumberGenerator& __urng,
1112  const param_type& __p)
1113  {
1114  result_type __q2 = __p.q() * __p.q();
1115  result_type __num = result_type(0.5L) * (result_type(1) + __q2);
1116  typename __gnu_cxx::arcsine_distribution<result_type>::param_type
1117  __pa(__num, __num / __q2);
1118  result_type __x = this->_M_ad(__pa, __urng);
1119  result_type __y = this->_M_ed(__urng);
1120  return (result_type(2) * __p.q() / (result_type(1) + __q2))
1121  * std::sqrt(__p.omega() * __x * __y);
1122  }
1123 
1124  template<typename _RealType>
1125  template<typename _OutputIterator,
1126  typename _UniformRandomNumberGenerator>
1127  void
1128  hoyt_distribution<_RealType>::
1129  __generate_impl(_OutputIterator __f, _OutputIterator __t,
1130  _UniformRandomNumberGenerator& __urng,
1131  const param_type& __p)
1132  {
1133  __glibcxx_function_requires(_OutputIteratorConcept<_OutputIterator,
1134  result_type>)
1135 
1136  result_type __2q = result_type(2) * __p.q();
1137  result_type __q2 = __p.q() * __p.q();
1138  result_type __q2p1 = result_type(1) + __q2;
1139  result_type __num = result_type(0.5L) * __q2p1;
1140  result_type __omega = __p.omega();
1141  typename __gnu_cxx::arcsine_distribution<result_type>::param_type
1142  __pa(__num, __num / __q2);
1143  while (__f != __t)
1144  {
1145  result_type __x = this->_M_ad(__pa, __urng);
1146  result_type __y = this->_M_ed(__urng);
1147  *__f++ = (__2q / __q2p1) * std::sqrt(__omega * __x * __y);
1148  }
1149  }
1150 
1151  template<typename _RealType, typename _CharT, typename _Traits>
1154  const hoyt_distribution<_RealType>& __x)
1155  {
1156  typedef std::basic_ostream<_CharT, _Traits> __ostream_type;
1157  typedef typename __ostream_type::ios_base __ios_base;
1158 
1159  const typename __ios_base::fmtflags __flags = __os.flags();
1160  const _CharT __fill = __os.fill();
1161  const std::streamsize __precision = __os.precision();
1162  const _CharT __space = __os.widen(' ');
1164  __os.fill(__space);
1166 
1167  __os << __x.q() << __space << __x.omega();
1168  __os << __space << __x._M_ad;
1169  __os << __space << __x._M_ed;
1170 
1171  __os.flags(__flags);
1172  __os.fill(__fill);
1173  __os.precision(__precision);
1174  return __os;
1175  }
1176 
1177  template<typename _RealType, typename _CharT, typename _Traits>
1180  hoyt_distribution<_RealType>& __x)
1181  {
1182  typedef std::basic_istream<_CharT, _Traits> __istream_type;
1183  typedef typename __istream_type::ios_base __ios_base;
1184 
1185  const typename __ios_base::fmtflags __flags = __is.flags();
1187 
1188  _RealType __q, __omega;
1189  __is >> __q >> __omega;
1190  __is >> __x._M_ad;
1191  __is >> __x._M_ed;
1192  __x.param(typename hoyt_distribution<_RealType>::
1193  param_type(__q, __omega));
1194 
1195  __is.flags(__flags);
1196  return __is;
1197  }
1198 
1199 
1200  template<typename _RealType>
1201  template<typename _OutputIterator,
1202  typename _UniformRandomNumberGenerator>
1203  void
1204  triangular_distribution<_RealType>::
1205  __generate_impl(_OutputIterator __f, _OutputIterator __t,
1206  _UniformRandomNumberGenerator& __urng,
1207  const param_type& __param)
1208  {
1209  __glibcxx_function_requires(_OutputIteratorConcept<_OutputIterator,
1210  result_type>)
1211 
1212  while (__f != __t)
1213  *__f++ = this->operator()(__urng, __param);
1214  }
1215 
1216  template<typename _RealType, typename _CharT, typename _Traits>
1217  std::basic_ostream<_CharT, _Traits>&
1218  operator<<(std::basic_ostream<_CharT, _Traits>& __os,
1219  const __gnu_cxx::triangular_distribution<_RealType>& __x)
1220  {
1221  typedef std::basic_ostream<_CharT, _Traits> __ostream_type;
1222  typedef typename __ostream_type::ios_base __ios_base;
1223 
1224  const typename __ios_base::fmtflags __flags = __os.flags();
1225  const _CharT __fill = __os.fill();
1226  const std::streamsize __precision = __os.precision();
1227  const _CharT __space = __os.widen(' ');
1229  __os.fill(__space);
1231 
1232  __os << __x.a() << __space << __x.b() << __space << __x.c();
1233 
1234  __os.flags(__flags);
1235  __os.fill(__fill);
1236  __os.precision(__precision);
1237  return __os;
1238  }
1239 
1240  template<typename _RealType, typename _CharT, typename _Traits>
1243  __gnu_cxx::triangular_distribution<_RealType>& __x)
1244  {
1245  typedef std::basic_istream<_CharT, _Traits> __istream_type;
1246  typedef typename __istream_type::ios_base __ios_base;
1247 
1248  const typename __ios_base::fmtflags __flags = __is.flags();
1250 
1251  _RealType __a, __b, __c;
1252  __is >> __a >> __b >> __c;
1253  __x.param(typename __gnu_cxx::triangular_distribution<_RealType>::
1254  param_type(__a, __b, __c));
1255 
1256  __is.flags(__flags);
1257  return __is;
1258  }
1259 
1260 
1261  template<typename _RealType>
1262  template<typename _UniformRandomNumberGenerator>
1263  typename von_mises_distribution<_RealType>::result_type
1264  von_mises_distribution<_RealType>::
1265  operator()(_UniformRandomNumberGenerator& __urng,
1266  const param_type& __p)
1267  {
1268  const result_type __pi
1269  = __gnu_cxx::__math_constants<result_type>::__pi;
1270  std::__detail::_Adaptor<_UniformRandomNumberGenerator, result_type>
1271  __aurng(__urng);
1272 
1273  result_type __f;
1274  while (1)
1275  {
1276  result_type __rnd = std::cos(__pi * __aurng());
1277  __f = (result_type(1) + __p._M_r * __rnd) / (__p._M_r + __rnd);
1278  result_type __c = __p._M_kappa * (__p._M_r - __f);
1279 
1280  result_type __rnd2 = __aurng();
1281  if (__c * (result_type(2) - __c) > __rnd2)
1282  break;
1283  if (std::log(__c / __rnd2) >= __c - result_type(1))
1284  break;
1285  }
1286 
1287  result_type __res = std::acos(__f);
1288 #if _GLIBCXX_USE_C99_MATH_TR1
1289  __res = std::copysign(__res, __aurng() - result_type(0.5));
1290 #else
1291  if (__aurng() < result_type(0.5))
1292  __res = -__res;
1293 #endif
1294  __res += __p._M_mu;
1295  if (__res > __pi)
1296  __res -= result_type(2) * __pi;
1297  else if (__res < -__pi)
1298  __res += result_type(2) * __pi;
1299  return __res;
1300  }
1301 
1302  template<typename _RealType>
1303  template<typename _OutputIterator,
1304  typename _UniformRandomNumberGenerator>
1305  void
1306  von_mises_distribution<_RealType>::
1307  __generate_impl(_OutputIterator __f, _OutputIterator __t,
1308  _UniformRandomNumberGenerator& __urng,
1309  const param_type& __param)
1310  {
1311  __glibcxx_function_requires(_OutputIteratorConcept<_OutputIterator,
1312  result_type>)
1313 
1314  while (__f != __t)
1315  *__f++ = this->operator()(__urng, __param);
1316  }
1317 
1318  template<typename _RealType, typename _CharT, typename _Traits>
1319  std::basic_ostream<_CharT, _Traits>&
1320  operator<<(std::basic_ostream<_CharT, _Traits>& __os,
1321  const __gnu_cxx::von_mises_distribution<_RealType>& __x)
1322  {
1323  typedef std::basic_ostream<_CharT, _Traits> __ostream_type;
1324  typedef typename __ostream_type::ios_base __ios_base;
1325 
1326  const typename __ios_base::fmtflags __flags = __os.flags();
1327  const _CharT __fill = __os.fill();
1328  const std::streamsize __precision = __os.precision();
1329  const _CharT __space = __os.widen(' ');
1331  __os.fill(__space);
1333 
1334  __os << __x.mu() << __space << __x.kappa();
1335 
1336  __os.flags(__flags);
1337  __os.fill(__fill);
1338  __os.precision(__precision);
1339  return __os;
1340  }
1341 
1342  template<typename _RealType, typename _CharT, typename _Traits>
1345  __gnu_cxx::von_mises_distribution<_RealType>& __x)
1346  {
1347  typedef std::basic_istream<_CharT, _Traits> __istream_type;
1348  typedef typename __istream_type::ios_base __ios_base;
1349 
1350  const typename __ios_base::fmtflags __flags = __is.flags();
1352 
1353  _RealType __mu, __kappa;
1354  __is >> __mu >> __kappa;
1355  __x.param(typename __gnu_cxx::von_mises_distribution<_RealType>::
1356  param_type(__mu, __kappa));
1357 
1358  __is.flags(__flags);
1359  return __is;
1360  }
1361 
1362 
1363  template<typename _UIntType>
1364  template<typename _UniformRandomNumberGenerator>
1365  typename hypergeometric_distribution<_UIntType>::result_type
1366  hypergeometric_distribution<_UIntType>::
1367  operator()(_UniformRandomNumberGenerator& __urng,
1368  const param_type& __param)
1369  {
1370  std::__detail::_Adaptor<_UniformRandomNumberGenerator, double>
1371  __aurng(__urng);
1372 
1373  result_type __a = __param.successful_size();
1374  result_type __b = __param.total_size();
1375  result_type __k = 0;
1376 
1377  if (__param.total_draws() < __param.total_size() / 2)
1378  {
1379  for (result_type __i = 0; __i < __param.total_draws(); ++__i)
1380  {
1381  if (__b * __aurng() < __a)
1382  {
1383  ++__k;
1384  if (__k == __param.successful_size())
1385  return __k;
1386  --__a;
1387  }
1388  --__b;
1389  }
1390  return __k;
1391  }
1392  else
1393  {
1394  for (result_type __i = 0; __i < __param.unsuccessful_size(); ++__i)
1395  {
1396  if (__b * __aurng() < __a)
1397  {
1398  ++__k;
1399  if (__k == __param.successful_size())
1400  return __param.successful_size() - __k;
1401  --__a;
1402  }
1403  --__b;
1404  }
1405  return __param.successful_size() - __k;
1406  }
1407  }
1408 
1409  template<typename _UIntType>
1410  template<typename _OutputIterator,
1411  typename _UniformRandomNumberGenerator>
1412  void
1413  hypergeometric_distribution<_UIntType>::
1414  __generate_impl(_OutputIterator __f, _OutputIterator __t,
1415  _UniformRandomNumberGenerator& __urng,
1416  const param_type& __param)
1417  {
1418  __glibcxx_function_requires(_OutputIteratorConcept<_OutputIterator,
1419  result_type>)
1420 
1421  while (__f != __t)
1422  *__f++ = this->operator()(__urng);
1423  }
1424 
1425  template<typename _UIntType, typename _CharT, typename _Traits>
1426  std::basic_ostream<_CharT, _Traits>&
1427  operator<<(std::basic_ostream<_CharT, _Traits>& __os,
1428  const __gnu_cxx::hypergeometric_distribution<_UIntType>& __x)
1429  {
1430  typedef std::basic_ostream<_CharT, _Traits> __ostream_type;
1431  typedef typename __ostream_type::ios_base __ios_base;
1432 
1433  const typename __ios_base::fmtflags __flags = __os.flags();
1434  const _CharT __fill = __os.fill();
1435  const std::streamsize __precision = __os.precision();
1436  const _CharT __space = __os.widen(' ');
1438  __os.fill(__space);
1440 
1441  __os << __x.total_size() << __space << __x.successful_size() << __space
1442  << __x.total_draws();
1443 
1444  __os.flags(__flags);
1445  __os.fill(__fill);
1446  __os.precision(__precision);
1447  return __os;
1448  }
1449 
1450  template<typename _UIntType, typename _CharT, typename _Traits>
1453  __gnu_cxx::hypergeometric_distribution<_UIntType>& __x)
1454  {
1455  typedef std::basic_istream<_CharT, _Traits> __istream_type;
1456  typedef typename __istream_type::ios_base __ios_base;
1457 
1458  const typename __ios_base::fmtflags __flags = __is.flags();
1460 
1461  _UIntType __total_size, __successful_size, __total_draws;
1462  __is >> __total_size >> __successful_size >> __total_draws;
1463  __x.param(typename __gnu_cxx::hypergeometric_distribution<_UIntType>::
1464  param_type(__total_size, __successful_size, __total_draws));
1465 
1466  __is.flags(__flags);
1467  return __is;
1468  }
1469 
1470 
1471  template<typename _RealType>
1472  template<typename _UniformRandomNumberGenerator>
1473  typename logistic_distribution<_RealType>::result_type
1474  logistic_distribution<_RealType>::
1475  operator()(_UniformRandomNumberGenerator& __urng,
1476  const param_type& __p)
1477  {
1478  std::__detail::_Adaptor<_UniformRandomNumberGenerator, result_type>
1479  __aurng(__urng);
1480 
1481  result_type __arg = result_type(1);
1482  while (__arg == result_type(1) || __arg == result_type(0))
1483  __arg = __aurng();
1484  return __p.a()
1485  + __p.b() * std::log(__arg / (result_type(1) - __arg));
1486  }
1487 
1488  template<typename _RealType>
1489  template<typename _OutputIterator,
1490  typename _UniformRandomNumberGenerator>
1491  void
1492  logistic_distribution<_RealType>::
1493  __generate_impl(_OutputIterator __f, _OutputIterator __t,
1494  _UniformRandomNumberGenerator& __urng,
1495  const param_type& __p)
1496  {
1497  __glibcxx_function_requires(_OutputIteratorConcept<_OutputIterator,
1498  result_type>)
1499 
1500  std::__detail::_Adaptor<_UniformRandomNumberGenerator, result_type>
1501  __aurng(__urng);
1502 
1503  while (__f != __t)
1504  {
1505  result_type __arg = result_type(1);
1506  while (__arg == result_type(1) || __arg == result_type(0))
1507  __arg = __aurng();
1508  *__f++ = __p.a()
1509  + __p.b() * std::log(__arg / (result_type(1) - __arg));
1510  }
1511  }
1512 
1513  template<typename _RealType, typename _CharT, typename _Traits>
1516  const logistic_distribution<_RealType>& __x)
1517  {
1518  typedef std::basic_ostream<_CharT, _Traits> __ostream_type;
1519  typedef typename __ostream_type::ios_base __ios_base;
1520 
1521  const typename __ios_base::fmtflags __flags = __os.flags();
1522  const _CharT __fill = __os.fill();
1523  const std::streamsize __precision = __os.precision();
1524  const _CharT __space = __os.widen(' ');
1526  __os.fill(__space);
1528 
1529  __os << __x.a() << __space << __x.b();
1530 
1531  __os.flags(__flags);
1532  __os.fill(__fill);
1533  __os.precision(__precision);
1534  return __os;
1535  }
1536 
1537  template<typename _RealType, typename _CharT, typename _Traits>
1540  logistic_distribution<_RealType>& __x)
1541  {
1542  typedef std::basic_istream<_CharT, _Traits> __istream_type;
1543  typedef typename __istream_type::ios_base __ios_base;
1544 
1545  const typename __ios_base::fmtflags __flags = __is.flags();
1547 
1548  _RealType __a, __b;
1549  __is >> __a >> __b;
1550  __x.param(typename logistic_distribution<_RealType>::
1551  param_type(__a, __b));
1552 
1553  __is.flags(__flags);
1554  return __is;
1555  }
1556 
1557 
1558  namespace {
1559 
1560  // Helper class for the uniform_on_sphere_distribution generation
1561  // function.
1562  template<std::size_t _Dimen, typename _RealType>
1563  class uniform_on_sphere_helper
1564  {
1565  typedef typename uniform_on_sphere_distribution<_Dimen, _RealType>::
1566  result_type result_type;
1567 
1568  public:
1569  template<typename _NormalDistribution,
1570  typename _UniformRandomNumberGenerator>
1571  result_type operator()(_NormalDistribution& __nd,
1572  _UniformRandomNumberGenerator& __urng)
1573  {
1574  result_type __ret;
1575  typename result_type::value_type __norm;
1576 
1577  do
1578  {
1579  auto __sum = _RealType(0);
1580 
1581  std::generate(__ret.begin(), __ret.end(),
1582  [&__nd, &__urng, &__sum](){
1583  _RealType __t = __nd(__urng);
1584  __sum += __t * __t;
1585  return __t; });
1586  __norm = std::sqrt(__sum);
1587  }
1588  while (__norm == _RealType(0) || ! __builtin_isfinite(__norm));
1589 
1590  std::transform(__ret.begin(), __ret.end(), __ret.begin(),
1591  [__norm](_RealType __val){ return __val / __norm; });
1592 
1593  return __ret;
1594  }
1595  };
1596 
1597 
1598  template<typename _RealType>
1599  class uniform_on_sphere_helper<2, _RealType>
1600  {
1601  typedef typename uniform_on_sphere_distribution<2, _RealType>::
1602  result_type result_type;
1603 
1604  public:
1605  template<typename _NormalDistribution,
1606  typename _UniformRandomNumberGenerator>
1607  result_type operator()(_NormalDistribution&,
1608  _UniformRandomNumberGenerator& __urng)
1609  {
1610  result_type __ret;
1611  _RealType __sq;
1612  std::__detail::_Adaptor<_UniformRandomNumberGenerator,
1613  _RealType> __aurng(__urng);
1614 
1615  do
1616  {
1617  __ret[0] = _RealType(2) * __aurng() - _RealType(1);
1618  __ret[1] = _RealType(2) * __aurng() - _RealType(1);
1619 
1620  __sq = __ret[0] * __ret[0] + __ret[1] * __ret[1];
1621  }
1622  while (__sq == _RealType(0) || __sq > _RealType(1));
1623 
1624 #if _GLIBCXX_USE_C99_MATH_TR1
1625  // Yes, we do not just use sqrt(__sq) because hypot() is more
1626  // accurate.
1627  auto __norm = std::hypot(__ret[0], __ret[1]);
1628 #else
1629  auto __norm = std::sqrt(__sq);
1630 #endif
1631  __ret[0] /= __norm;
1632  __ret[1] /= __norm;
1633 
1634  return __ret;
1635  }
1636  };
1637 
1638  }
1639 
1640 
1641  template<std::size_t _Dimen, typename _RealType>
1642  template<typename _UniformRandomNumberGenerator>
1643  typename uniform_on_sphere_distribution<_Dimen, _RealType>::result_type
1644  uniform_on_sphere_distribution<_Dimen, _RealType>::
1645  operator()(_UniformRandomNumberGenerator& __urng,
1646  const param_type& __p)
1647  {
1648  uniform_on_sphere_helper<_Dimen, _RealType> __helper;
1649  return __helper(_M_nd, __urng);
1650  }
1651 
1652  template<std::size_t _Dimen, typename _RealType>
1653  template<typename _OutputIterator,
1654  typename _UniformRandomNumberGenerator>
1655  void
1656  uniform_on_sphere_distribution<_Dimen, _RealType>::
1657  __generate_impl(_OutputIterator __f, _OutputIterator __t,
1658  _UniformRandomNumberGenerator& __urng,
1659  const param_type& __param)
1660  {
1661  __glibcxx_function_requires(_OutputIteratorConcept<_OutputIterator,
1662  result_type>)
1663 
1664  while (__f != __t)
1665  *__f++ = this->operator()(__urng, __param);
1666  }
1667 
1668  template<std::size_t _Dimen, typename _RealType, typename _CharT,
1669  typename _Traits>
1670  std::basic_ostream<_CharT, _Traits>&
1671  operator<<(std::basic_ostream<_CharT, _Traits>& __os,
1672  const __gnu_cxx::uniform_on_sphere_distribution<_Dimen,
1673  _RealType>& __x)
1674  {
1675  return __os << __x._M_nd;
1676  }
1677 
1678  template<std::size_t _Dimen, typename _RealType, typename _CharT,
1679  typename _Traits>
1682  __gnu_cxx::uniform_on_sphere_distribution<_Dimen,
1683  _RealType>& __x)
1684  {
1685  return __is >> __x._M_nd;
1686  }
1687 
1688 
1689  namespace {
1690 
1691  // Helper class for the uniform_inside_sphere_distribution generation
1692  // function.
1693  template<std::size_t _Dimen, bool _SmallDimen, typename _RealType>
1694  class uniform_inside_sphere_helper;
1695 
1696  template<std::size_t _Dimen, typename _RealType>
1697  class uniform_inside_sphere_helper<_Dimen, false, _RealType>
1698  {
1699  using result_type
1700  = typename uniform_inside_sphere_distribution<_Dimen, _RealType>::
1701  result_type;
1702 
1703  public:
1704  template<typename _UniformOnSphereDistribution,
1705  typename _UniformRandomNumberGenerator>
1706  result_type
1707  operator()(_UniformOnSphereDistribution& __uosd,
1708  _UniformRandomNumberGenerator& __urng,
1709  _RealType __radius)
1710  {
1711  std::__detail::_Adaptor<_UniformRandomNumberGenerator,
1712  _RealType> __aurng(__urng);
1713 
1714  _RealType __pow = 1 / _RealType(_Dimen);
1715  _RealType __urt = __radius * std::pow(__aurng(), __pow);
1716  result_type __ret = __uosd(__aurng);
1717 
1718  std::transform(__ret.begin(), __ret.end(), __ret.begin(),
1719  [__urt](_RealType __val)
1720  { return __val * __urt; });
1721 
1722  return __ret;
1723  }
1724  };
1725 
1726  // Helper class for the uniform_inside_sphere_distribution generation
1727  // function specialized for small dimensions.
1728  template<std::size_t _Dimen, typename _RealType>
1729  class uniform_inside_sphere_helper<_Dimen, true, _RealType>
1730  {
1731  using result_type
1732  = typename uniform_inside_sphere_distribution<_Dimen, _RealType>::
1733  result_type;
1734 
1735  public:
1736  template<typename _UniformOnSphereDistribution,
1737  typename _UniformRandomNumberGenerator>
1738  result_type
1739  operator()(_UniformOnSphereDistribution&,
1740  _UniformRandomNumberGenerator& __urng,
1741  _RealType __radius)
1742  {
1743  result_type __ret;
1744  _RealType __sq;
1745  _RealType __radsq = __radius * __radius;
1746  std::__detail::_Adaptor<_UniformRandomNumberGenerator,
1747  _RealType> __aurng(__urng);
1748 
1749  do
1750  {
1751  __sq = _RealType(0);
1752  for (int i = 0; i < _Dimen; ++i)
1753  {
1754  __ret[i] = _RealType(2) * __aurng() - _RealType(1);
1755  __sq += __ret[i] * __ret[i];
1756  }
1757  }
1758  while (__sq > _RealType(1));
1759 
1760  for (int i = 0; i < _Dimen; ++i)
1761  __ret[i] *= __radius;
1762 
1763  return __ret;
1764  }
1765  };
1766  } // namespace
1767 
1768  //
1769  // Experiments have shown that rejection is more efficient than transform
1770  // for dimensions less than 8.
1771  //
1772  template<std::size_t _Dimen, typename _RealType>
1773  template<typename _UniformRandomNumberGenerator>
1774  typename uniform_inside_sphere_distribution<_Dimen, _RealType>::result_type
1775  uniform_inside_sphere_distribution<_Dimen, _RealType>::
1776  operator()(_UniformRandomNumberGenerator& __urng,
1777  const param_type& __p)
1778  {
1779  uniform_inside_sphere_helper<_Dimen, _Dimen < 8, _RealType> __helper;
1780  return __helper(_M_uosd, __urng, __p.radius());
1781  }
1782 
1783  template<std::size_t _Dimen, typename _RealType>
1784  template<typename _OutputIterator,
1785  typename _UniformRandomNumberGenerator>
1786  void
1787  uniform_inside_sphere_distribution<_Dimen, _RealType>::
1788  __generate_impl(_OutputIterator __f, _OutputIterator __t,
1789  _UniformRandomNumberGenerator& __urng,
1790  const param_type& __param)
1791  {
1792  __glibcxx_function_requires(_OutputIteratorConcept<_OutputIterator,
1793  result_type>)
1794 
1795  while (__f != __t)
1796  *__f++ = this->operator()(__urng, __param);
1797  }
1798 
1799  template<std::size_t _Dimen, typename _RealType, typename _CharT,
1800  typename _Traits>
1801  std::basic_ostream<_CharT, _Traits>&
1802  operator<<(std::basic_ostream<_CharT, _Traits>& __os,
1803  const __gnu_cxx::uniform_inside_sphere_distribution<_Dimen,
1804  _RealType>& __x)
1805  {
1806  typedef std::basic_ostream<_CharT, _Traits> __ostream_type;
1807  typedef typename __ostream_type::ios_base __ios_base;
1808 
1809  const typename __ios_base::fmtflags __flags = __os.flags();
1810  const _CharT __fill = __os.fill();
1811  const std::streamsize __precision = __os.precision();
1812  const _CharT __space = __os.widen(' ');
1814  __os.fill(__space);
1816 
1817  __os << __x.radius() << __space << __x._M_uosd;
1818 
1819  __os.flags(__flags);
1820  __os.fill(__fill);
1821  __os.precision(__precision);
1822 
1823  return __os;
1824  }
1825 
1826  template<std::size_t _Dimen, typename _RealType, typename _CharT,
1827  typename _Traits>
1830  __gnu_cxx::uniform_inside_sphere_distribution<_Dimen,
1831  _RealType>& __x)
1832  {
1833  typedef std::basic_istream<_CharT, _Traits> __istream_type;
1834  typedef typename __istream_type::ios_base __ios_base;
1835 
1836  const typename __ios_base::fmtflags __flags = __is.flags();
1838 
1839  _RealType __radius_val;
1840  __is >> __radius_val >> __x._M_uosd;
1841  __x.param(typename uniform_inside_sphere_distribution<_Dimen, _RealType>::
1842  param_type(__radius_val));
1843 
1844  __is.flags(__flags);
1845 
1846  return __is;
1847  }
1848 
1849 _GLIBCXX_END_NAMESPACE_VERSION
1850 } // namespace __gnu_cxx
1851 
1852 
1853 #endif // _EXT_RANDOM_TCC
complex< _Tp > sin(const complex< _Tp > &)
Return complex sine of z.
Definition: complex:859
complex< _Tp > log(const complex< _Tp > &)
Return complex natural logarithm of z.
Definition: complex:824
complex< _Tp > exp(const complex< _Tp > &)
Return complex base e exponential of z.
Definition: complex:797
complex< _Tp > pow(const complex< _Tp > &, int)
Return x to the y'th power.
Definition: complex:1019
complex< _Tp > cos(const complex< _Tp > &)
Return complex cosine of z.
Definition: complex:741
complex< _Tp > sqrt(const complex< _Tp > &)
Return complex square root of z.
Definition: complex:933
ISO C++ entities toplevel namespace is std.
std::basic_istream< _CharT, _Traits > & operator>>(std::basic_istream< _CharT, _Traits > &__is, bitset< _Nb > &__x)
Global I/O operators for bitsets.
Definition: bitset:1472
ptrdiff_t streamsize
Integral type for I/O operation counts and buffer sizes.
Definition: postypes.h:68
std::basic_ostream< _CharT, _Traits > & operator<<(std::basic_ostream< _CharT, _Traits > &__os, const bitset< _Nb > &__x)
Global I/O operators for bitsets.
Definition: bitset:1540
ios_base & scientific(ios_base &__base)
Calls base.setf(ios_base::scientific, ios_base::floatfield).
Definition: ios_base.h:1079
ios_base & dec(ios_base &__base)
Calls base.setf(ios_base::dec, ios_base::basefield).
Definition: ios_base.h:1046
ios_base & left(ios_base &__base)
Calls base.setf(ios_base::left, ios_base::adjustfield).
Definition: ios_base.h:1029
ios_base & skipws(ios_base &__base)
Calls base.setf(ios_base::skipws).
Definition: ios_base.h:972
std::complex< _Tp > acos(const std::complex< _Tp > &)
acos(__z) [8.1.2].
Definition: complex:1638
ios_base & fixed(ios_base &__base)
Calls base.setf(ios_base::fixed, ios_base::floatfield).
Definition: ios_base.h:1071
GNU extensions for public use.
A standard container for storing a fixed size sequence of elements.
Definition: array:96
char_type widen(char __c) const
Widens characters.
Definition: basic_ios.h:449
char_type fill() const
Retrieves the empty character.
Definition: basic_ios.h:370
Template class basic_istream.
Definition: istream:59
Template class basic_ostream.
Definition: ostream:59
Properties of fundamental types.
Definition: limits:313
streamsize precision() const
Flags access.
Definition: ios_base.h:719
fmtflags flags() const
Access to format flags.
Definition: ios_base.h:649