Clang Project

include/c++/7/bits/random.tcc
1// random number generation (out of line) -*- C++ -*-
2
3// Copyright (C) 2009-2017 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 bits/random.tcc
26 *  This is an internal header file, included by other library headers.
27 *  Do not attempt to use it directly. @headername{random}
28 */
29
30#ifndef _RANDOM_TCC
31#define _RANDOM_TCC 1
32
33#include <numeric> // std::accumulate and std::partial_sum
34
35namespace std _GLIBCXX_VISIBILITY(default)
36{
37  /*
38   * (Further) implementation-space details.
39   */
40  namespace __detail
41  {
42  _GLIBCXX_BEGIN_NAMESPACE_VERSION
43
44    // General case for x = (ax + c) mod m -- use Schrage's algorithm
45    // to avoid integer overflow.
46    //
47    // Preconditions:  a > 0, m > 0.
48    //
49    // Note: only works correctly for __m % __a < __m / __a.
50    template<typename _Tp, _Tp __m, _Tp __a, _Tp __c>
51      _Tp
52      _Mod<_Tp, __m, __a, __c, false, true>::
53      __calc(_Tp __x)
54      {
55 if (__a == 1)
56   __x %= __m;
57 else
58   {
59     static const _Tp __q = __m / __a;
60     static const _Tp __r = __m % __a;
61
62     _Tp __t1 = __a * (__x % __q);
63     _Tp __t2 = __r * (__x / __q);
64     if (__t1 >= __t2)
65       __x = __t1 - __t2;
66     else
67       __x = __m - __t2 + __t1;
68   }
69
70 if (__c != 0)
71   {
72     const _Tp __d = __m - __x;
73     if (__d > __c)
74       __x += __c;
75     else
76       __x = __c - __d;
77   }
78 return __x;
79      }
80
81    template<typename _InputIterator, typename _OutputIterator,
82      typename _Tp>
83      _OutputIterator
84      __normalize(_InputIterator __first, _InputIterator __last,
85   _OutputIterator __result, const _Tp& __factor)
86      {
87 for (; __first != __last; ++__first, ++__result)
88   *__result = *__first / __factor;
89 return __result;
90      }
91
92  _GLIBCXX_END_NAMESPACE_VERSION
93  } // namespace __detail
94
95_GLIBCXX_BEGIN_NAMESPACE_VERSION
96
97  template<typename _UIntType, _UIntType __a, _UIntType __c, _UIntType __m>
98    constexpr _UIntType
99    linear_congruential_engine<_UIntType, __a, __c, __m>::multiplier;
100
101  template<typename _UIntType, _UIntType __a, _UIntType __c, _UIntType __m>
102    constexpr _UIntType
103    linear_congruential_engine<_UIntType, __a, __c, __m>::increment;
104
105  template<typename _UIntType, _UIntType __a, _UIntType __c, _UIntType __m>
106    constexpr _UIntType
107    linear_congruential_engine<_UIntType, __a, __c, __m>::modulus;
108
109  template<typename _UIntType, _UIntType __a, _UIntType __c, _UIntType __m>
110    constexpr _UIntType
111    linear_congruential_engine<_UIntType, __a, __c, __m>::default_seed;
112
113  /**
114   * Seeds the LCR with integral value @p __s, adjusted so that the
115   * ring identity is never a member of the convergence set.
116   */
117  template<typename _UIntType, _UIntType __a, _UIntType __c, _UIntType __m>
118    void
119    linear_congruential_engine<_UIntType, __a, __c, __m>::
120    seed(result_type __s)
121    {
122      if ((__detail::__mod<_UIntType, __m>(__c) == 0)
123   && (__detail::__mod<_UIntType, __m>(__s) == 0))
124 _M_x = 1;
125      else
126 _M_x = __detail::__mod<_UIntType, __m>(__s);
127    }
128
129  /**
130   * Seeds the LCR engine with a value generated by @p __q.
131   */
132  template<typename _UIntType, _UIntType __a, _UIntType __c, _UIntType __m>
133    template<typename _Sseq>
134      typename std::enable_if<std::is_class<_Sseq>::value>::type
135      linear_congruential_engine<_UIntType, __a, __c, __m>::
136      seed(_Sseq& __q)
137      {
138 const _UIntType __k0 = __m == 0 ? std::numeric_limits<_UIntType>::digits
139                                 : std::__lg(__m);
140 const _UIntType __k = (__k0 + 31) / 32;
141 uint_least32_t __arr[__k + 3];
142 __q.generate(__arr + 0, __arr + __k + 3);
143 _UIntType __factor = 1u;
144 _UIntType __sum = 0u;
145 for (size_t __j = 0; __j < __k; ++__j)
146   {
147     __sum += __arr[__j + 3] * __factor;
148     __factor *= __detail::_Shift<_UIntType, 32>::__value;
149   }
150 seed(__sum);
151      }
152
153  template<typename _UIntType, _UIntType __a, _UIntType __c, _UIntType __m,
154    typename _CharT, typename _Traits>
155    std::basic_ostream<_CharT, _Traits>&
156    operator<<(std::basic_ostream<_CharT, _Traits>& __os,
157        const linear_congruential_engine<_UIntType,
158 __a, __c, __m>& __lcr)
159    {
160      typedef std::basic_ostream<_CharT, _Traits>  __ostream_type;
161      typedef typename __ostream_type::ios_base    __ios_base;
162
163      const typename __ios_base::fmtflags __flags = __os.flags();
164      const _CharT __fill = __os.fill();
165      __os.flags(__ios_base::dec | __ios_base::fixed | __ios_base::left);
166      __os.fill(__os.widen(' '));
167
168      __os << __lcr._M_x;
169
170      __os.flags(__flags);
171      __os.fill(__fill);
172      return __os;
173    }
174
175  template<typename _UIntType, _UIntType __a, _UIntType __c, _UIntType __m,
176    typename _CharT, typename _Traits>
177    std::basic_istream<_CharT, _Traits>&
178    operator>>(std::basic_istream<_CharT, _Traits>& __is,
179        linear_congruential_engine<_UIntType, __a, __c, __m>& __lcr)
180    {
181      typedef std::basic_istream<_CharT, _Traits>  __istream_type;
182      typedef typename __istream_type::ios_base    __ios_base;
183
184      const typename __ios_base::fmtflags __flags = __is.flags();
185      __is.flags(__ios_base::dec);
186
187      __is >> __lcr._M_x;
188
189      __is.flags(__flags);
190      return __is;
191    }
192
193
194  template<typename _UIntType,
195    size_t __w, size_t __n, size_t __m, size_t __r,
196    _UIntType __a, size_t __u, _UIntType __d, size_t __s,
197    _UIntType __b, size_t __t, _UIntType __c, size_t __l,
198    _UIntType __f>
199    constexpr size_t
200    mersenne_twister_engine<_UIntType, __w, __n, __m, __r, __a, __u, __d,
201     __s, __b, __t, __c, __l, __f>::word_size;
202
203  template<typename _UIntType,
204    size_t __w, size_t __n, size_t __m, size_t __r,
205    _UIntType __a, size_t __u, _UIntType __d, size_t __s,
206    _UIntType __b, size_t __t, _UIntType __c, size_t __l,
207    _UIntType __f>
208    constexpr size_t
209    mersenne_twister_engine<_UIntType, __w, __n, __m, __r, __a, __u, __d,
210     __s, __b, __t, __c, __l, __f>::state_size;
211
212  template<typename _UIntType,
213    size_t __w, size_t __n, size_t __m, size_t __r,
214    _UIntType __a, size_t __u, _UIntType __d, size_t __s,
215    _UIntType __b, size_t __t, _UIntType __c, size_t __l,
216    _UIntType __f>
217    constexpr size_t
218    mersenne_twister_engine<_UIntType, __w, __n, __m, __r, __a, __u, __d,
219     __s, __b, __t, __c, __l, __f>::shift_size;
220
221  template<typename _UIntType,
222    size_t __w, size_t __n, size_t __m, size_t __r,
223    _UIntType __a, size_t __u, _UIntType __d, size_t __s,
224    _UIntType __b, size_t __t, _UIntType __c, size_t __l,
225    _UIntType __f>
226    constexpr size_t
227    mersenne_twister_engine<_UIntType, __w, __n, __m, __r, __a, __u, __d,
228     __s, __b, __t, __c, __l, __f>::mask_bits;
229
230  template<typename _UIntType,
231    size_t __w, size_t __n, size_t __m, size_t __r,
232    _UIntType __a, size_t __u, _UIntType __d, size_t __s,
233    _UIntType __b, size_t __t, _UIntType __c, size_t __l,
234    _UIntType __f>
235    constexpr _UIntType
236    mersenne_twister_engine<_UIntType, __w, __n, __m, __r, __a, __u, __d,
237     __s, __b, __t, __c, __l, __f>::xor_mask;
238
239  template<typename _UIntType,
240    size_t __w, size_t __n, size_t __m, size_t __r,
241    _UIntType __a, size_t __u, _UIntType __d, size_t __s,
242    _UIntType __b, size_t __t, _UIntType __c, size_t __l,
243    _UIntType __f>
244    constexpr size_t
245    mersenne_twister_engine<_UIntType, __w, __n, __m, __r, __a, __u, __d,
246     __s, __b, __t, __c, __l, __f>::tempering_u;
247   
248  template<typename _UIntType,
249    size_t __w, size_t __n, size_t __m, size_t __r,
250    _UIntType __a, size_t __u, _UIntType __d, size_t __s,
251    _UIntType __b, size_t __t, _UIntType __c, size_t __l,
252    _UIntType __f>
253    constexpr _UIntType
254    mersenne_twister_engine<_UIntType, __w, __n, __m, __r, __a, __u, __d,
255     __s, __b, __t, __c, __l, __f>::tempering_d;
256
257  template<typename _UIntType,
258    size_t __w, size_t __n, size_t __m, size_t __r,
259    _UIntType __a, size_t __u, _UIntType __d, size_t __s,
260    _UIntType __b, size_t __t, _UIntType __c, size_t __l,
261    _UIntType __f>
262    constexpr size_t
263    mersenne_twister_engine<_UIntType, __w, __n, __m, __r, __a, __u, __d,
264     __s, __b, __t, __c, __l, __f>::tempering_s;
265
266  template<typename _UIntType,
267    size_t __w, size_t __n, size_t __m, size_t __r,
268    _UIntType __a, size_t __u, _UIntType __d, size_t __s,
269    _UIntType __b, size_t __t, _UIntType __c, size_t __l,
270    _UIntType __f>
271    constexpr _UIntType
272    mersenne_twister_engine<_UIntType, __w, __n, __m, __r, __a, __u, __d,
273     __s, __b, __t, __c, __l, __f>::tempering_b;
274
275  template<typename _UIntType,
276    size_t __w, size_t __n, size_t __m, size_t __r,
277    _UIntType __a, size_t __u, _UIntType __d, size_t __s,
278    _UIntType __b, size_t __t, _UIntType __c, size_t __l,
279    _UIntType __f>
280    constexpr size_t
281    mersenne_twister_engine<_UIntType, __w, __n, __m, __r, __a, __u, __d,
282     __s, __b, __t, __c, __l, __f>::tempering_t;
283
284  template<typename _UIntType,
285    size_t __w, size_t __n, size_t __m, size_t __r,
286    _UIntType __a, size_t __u, _UIntType __d, size_t __s,
287    _UIntType __b, size_t __t, _UIntType __c, size_t __l,
288    _UIntType __f>
289    constexpr _UIntType
290    mersenne_twister_engine<_UIntType, __w, __n, __m, __r, __a, __u, __d,
291     __s, __b, __t, __c, __l, __f>::tempering_c;
292
293  template<typename _UIntType,
294    size_t __w, size_t __n, size_t __m, size_t __r,
295    _UIntType __a, size_t __u, _UIntType __d, size_t __s,
296    _UIntType __b, size_t __t, _UIntType __c, size_t __l,
297    _UIntType __f>
298    constexpr size_t
299    mersenne_twister_engine<_UIntType, __w, __n, __m, __r, __a, __u, __d,
300     __s, __b, __t, __c, __l, __f>::tempering_l;
301
302  template<typename _UIntType,
303    size_t __w, size_t __n, size_t __m, size_t __r,
304    _UIntType __a, size_t __u, _UIntType __d, size_t __s,
305    _UIntType __b, size_t __t, _UIntType __c, size_t __l,
306    _UIntType __f>
307    constexpr _UIntType
308    mersenne_twister_engine<_UIntType, __w, __n, __m, __r, __a, __u, __d,
309     __s, __b, __t, __c, __l, __f>::
310                                              initialization_multiplier;
311
312  template<typename _UIntType,
313    size_t __w, size_t __n, size_t __m, size_t __r,
314    _UIntType __a, size_t __u, _UIntType __d, size_t __s,
315    _UIntType __b, size_t __t, _UIntType __c, size_t __l,
316    _UIntType __f>
317    constexpr _UIntType
318    mersenne_twister_engine<_UIntType, __w, __n, __m, __r, __a, __u, __d,
319     __s, __b, __t, __c, __l, __f>::default_seed;
320
321  template<typename _UIntType,
322    size_t __w, size_t __n, size_t __m, size_t __r,
323    _UIntType __a, size_t __u, _UIntType __d, size_t __s,
324    _UIntType __b, size_t __t, _UIntType __c, size_t __l,
325    _UIntType __f>
326    void
327    mersenne_twister_engine<_UIntType, __w, __n, __m, __r, __a, __u, __d,
328     __s, __b, __t, __c, __l, __f>::
329    seed(result_type __sd)
330    {
331      _M_x[0] = __detail::__mod<_UIntType,
332 __detail::_Shift<_UIntType, __w>::__value>(__sd);
333
334      for (size_t __i = 1; __i < state_size; ++__i)
335 {
336   _UIntType __x = _M_x[__i - 1];
337   __x ^= __x >> (__w - 2);
338   __x *= __f;
339   __x += __detail::__mod<_UIntType, __n>(__i);
340   _M_x[__i] = __detail::__mod<_UIntType,
341     __detail::_Shift<_UIntType, __w>::__value>(__x);
342 }
343      _M_p = state_size;
344    }
345
346  template<typename _UIntType,
347    size_t __w, size_t __n, size_t __m, size_t __r,
348    _UIntType __a, size_t __u, _UIntType __d, size_t __s,
349    _UIntType __b, size_t __t, _UIntType __c, size_t __l,
350    _UIntType __f>
351    template<typename _Sseq>
352      typename std::enable_if<std::is_class<_Sseq>::value>::type
353      mersenne_twister_engine<_UIntType, __w, __n, __m, __r, __a, __u, __d,
354       __s, __b, __t, __c, __l, __f>::
355      seed(_Sseq& __q)
356      {
357 const _UIntType __upper_mask = (~_UIntType()) << __r;
358 const size_t __k = (__w + 31) / 32;
359 uint_least32_t __arr[__n * __k];
360 __q.generate(__arr + 0, __arr + __n * __k);
361
362 bool __zero = true;
363 for (size_t __i = 0; __i < state_size; ++__i)
364   {
365     _UIntType __factor = 1u;
366     _UIntType __sum = 0u;
367     for (size_t __j = 0; __j < __k; ++__j)
368       {
369 __sum += __arr[__k * __i + __j] * __factor;
370 __factor *= __detail::_Shift<_UIntType, 32>::__value;
371       }
372     _M_x[__i] = __detail::__mod<_UIntType,
373       __detail::_Shift<_UIntType, __w>::__value>(__sum);
374
375     if (__zero)
376       {
377 if (__i == 0)
378   {
379     if ((_M_x[0] & __upper_mask) != 0u)
380       __zero = false;
381   }
382 else if (_M_x[__i] != 0u)
383   __zero = false;
384       }
385   }
386        if (__zero)
387          _M_x[0] = __detail::_Shift<_UIntType, __w - 1>::__value;
388 _M_p = state_size;
389      }
390
391  template<typename _UIntType, size_t __w,
392    size_t __n, size_t __m, size_t __r,
393    _UIntType __a, size_t __u, _UIntType __d, size_t __s,
394    _UIntType __b, size_t __t, _UIntType __c, size_t __l,
395    _UIntType __f>
396    void
397    mersenne_twister_engine<_UIntType, __w, __n, __m, __r, __a, __u, __d,
398     __s, __b, __t, __c, __l, __f>::
399    _M_gen_rand(void)
400    {
401      const _UIntType __upper_mask = (~_UIntType()) << __r;
402      const _UIntType __lower_mask = ~__upper_mask;
403
404      for (size_t __k = 0; __k < (__n - __m); ++__k)
405        {
406   _UIntType __y = ((_M_x[__k] & __upper_mask)
407    | (_M_x[__k + 1] & __lower_mask));
408   _M_x[__k] = (_M_x[__k + __m] ^ (__y >> 1)
409        ^ ((__y & 0x01) ? __a : 0));
410        }
411
412      for (size_t __k = (__n - __m); __k < (__n - 1); ++__k)
413 {
414   _UIntType __y = ((_M_x[__k] & __upper_mask)
415    | (_M_x[__k + 1] & __lower_mask));
416   _M_x[__k] = (_M_x[__k + (__m - __n)] ^ (__y >> 1)
417        ^ ((__y & 0x01) ? __a : 0));
418 }
419
420      _UIntType __y = ((_M_x[__n - 1] & __upper_mask)
421        | (_M_x[0] & __lower_mask));
422      _M_x[__n - 1] = (_M_x[__m - 1] ^ (__y >> 1)
423        ^ ((__y & 0x01) ? __a : 0));
424      _M_p = 0;
425    }
426
427  template<typename _UIntType, size_t __w,
428    size_t __n, size_t __m, size_t __r,
429    _UIntType __a, size_t __u, _UIntType __d, size_t __s,
430    _UIntType __b, size_t __t, _UIntType __c, size_t __l,
431    _UIntType __f>
432    void
433    mersenne_twister_engine<_UIntType, __w, __n, __m, __r, __a, __u, __d,
434     __s, __b, __t, __c, __l, __f>::
435    discard(unsigned long long __z)
436    {
437      while (__z > state_size - _M_p)
438 {
439   __z -= state_size - _M_p;
440   _M_gen_rand();
441 }
442      _M_p += __z;
443    }
444
445  template<typename _UIntType, size_t __w,
446    size_t __n, size_t __m, size_t __r,
447    _UIntType __a, size_t __u, _UIntType __d, size_t __s,
448    _UIntType __b, size_t __t, _UIntType __c, size_t __l,
449    _UIntType __f>
450    typename
451    mersenne_twister_engine<_UIntType, __w, __n, __m, __r, __a, __u, __d,
452     __s, __b, __t, __c, __l, __f>::result_type
453    mersenne_twister_engine<_UIntType, __w, __n, __m, __r, __a, __u, __d,
454     __s, __b, __t, __c, __l, __f>::
455    operator()()
456    {
457      // Reload the vector - cost is O(n) amortized over n calls.
458      if (_M_p >= state_size)
459 _M_gen_rand();
460
461      // Calculate o(x(i)).
462      result_type __z = _M_x[_M_p++];
463      __z ^= (__z >> __u) & __d;
464      __z ^= (__z << __s) & __b;
465      __z ^= (__z << __t) & __c;
466      __z ^= (__z >> __l);
467
468      return __z;
469    }
470
471  template<typename _UIntType, size_t __w,
472    size_t __n, size_t __m, size_t __r,
473    _UIntType __a, size_t __u, _UIntType __d, size_t __s,
474    _UIntType __b, size_t __t, _UIntType __c, size_t __l,
475    _UIntType __f, typename _CharT, typename _Traits>
476    std::basic_ostream<_CharT, _Traits>&
477    operator<<(std::basic_ostream<_CharT, _Traits>& __os,
478        const mersenne_twister_engine<_UIntType, __w, __n, __m,
479        __r, __a, __u, __d, __s, __b, __t, __c, __l, __f>& __x)
480    {
481      typedef std::basic_ostream<_CharT, _Traits>  __ostream_type;
482      typedef typename __ostream_type::ios_base    __ios_base;
483
484      const typename __ios_base::fmtflags __flags = __os.flags();
485      const _CharT __fill = __os.fill();
486      const _CharT __space = __os.widen(' ');
487      __os.flags(__ios_base::dec | __ios_base::fixed | __ios_base::left);
488      __os.fill(__space);
489
490      for (size_t __i = 0; __i < __n; ++__i)
491 __os << __x._M_x[__i] << __space;
492      __os << __x._M_p;
493
494      __os.flags(__flags);
495      __os.fill(__fill);
496      return __os;
497    }
498
499  template<typename _UIntType, size_t __w,
500    size_t __n, size_t __m, size_t __r,
501    _UIntType __a, size_t __u, _UIntType __d, size_t __s,
502    _UIntType __b, size_t __t, _UIntType __c, size_t __l,
503    _UIntType __f, typename _CharT, typename _Traits>
504    std::basic_istream<_CharT, _Traits>&
505    operator>>(std::basic_istream<_CharT, _Traits>& __is,
506        mersenne_twister_engine<_UIntType, __w, __n, __m,
507        __r, __a, __u, __d, __s, __b, __t, __c, __l, __f>& __x)
508    {
509      typedef std::basic_istream<_CharT, _Traits>  __istream_type;
510      typedef typename __istream_type::ios_base    __ios_base;
511
512      const typename __ios_base::fmtflags __flags = __is.flags();
513      __is.flags(__ios_base::dec | __ios_base::skipws);
514
515      for (size_t __i = 0; __i < __n; ++__i)
516 __is >> __x._M_x[__i];
517      __is >> __x._M_p;
518
519      __is.flags(__flags);
520      return __is;
521    }
522
523
524  template<typename _UIntType, size_t __w, size_t __s, size_t __r>
525    constexpr size_t
526    subtract_with_carry_engine<_UIntType, __w, __s, __r>::word_size;
527
528  template<typename _UIntType, size_t __w, size_t __s, size_t __r>
529    constexpr size_t
530    subtract_with_carry_engine<_UIntType, __w, __s, __r>::short_lag;
531
532  template<typename _UIntType, size_t __w, size_t __s, size_t __r>
533    constexpr size_t
534    subtract_with_carry_engine<_UIntType, __w, __s, __r>::long_lag;
535
536  template<typename _UIntType, size_t __w, size_t __s, size_t __r>
537    constexpr _UIntType
538    subtract_with_carry_engine<_UIntType, __w, __s, __r>::default_seed;
539
540  template<typename _UIntType, size_t __w, size_t __s, size_t __r>
541    void
542    subtract_with_carry_engine<_UIntType, __w, __s, __r>::
543    seed(result_type __value)
544    {
545      std::linear_congruential_engine<result_type, 40014u, 0u, 2147483563u>
546 __lcg(__value == 0u ? default_seed : __value);
547
548      const size_t __n = (__w + 31) / 32;
549
550      for (size_t __i = 0; __i < long_lag; ++__i)
551 {
552   _UIntType __sum = 0u;
553   _UIntType __factor = 1u;
554   for (size_t __j = 0; __j < __n; ++__j)
555     {
556       __sum += __detail::__mod<uint_least32_t,
557        __detail::_Shift<uint_least32_t, 32>::__value>
558  (__lcg()) * __factor;
559       __factor *= __detail::_Shift<_UIntType, 32>::__value;
560     }
561   _M_x[__i] = __detail::__mod<_UIntType,
562     __detail::_Shift<_UIntType, __w>::__value>(__sum);
563 }
564      _M_carry = (_M_x[long_lag - 1] == 0) ? 1 : 0;
565      _M_p = 0;
566    }
567
568  template<typename _UIntType, size_t __w, size_t __s, size_t __r>
569    template<typename _Sseq>
570      typename std::enable_if<std::is_class<_Sseq>::value>::type
571      subtract_with_carry_engine<_UIntType, __w, __s, __r>::
572      seed(_Sseq& __q)
573      {
574 const size_t __k = (__w + 31) / 32;
575 uint_least32_t __arr[__r * __k];
576 __q.generate(__arr + 0, __arr + __r * __k);
577
578 for (size_t __i = 0; __i < long_lag; ++__i)
579   {
580     _UIntType __sum = 0u;
581     _UIntType __factor = 1u;
582     for (size_t __j = 0; __j < __k; ++__j)
583       {
584 __sum += __arr[__k * __i + __j] * __factor;
585 __factor *= __detail::_Shift<_UIntType, 32>::__value;
586       }
587     _M_x[__i] = __detail::__mod<_UIntType,
588       __detail::_Shift<_UIntType, __w>::__value>(__sum);
589   }
590 _M_carry = (_M_x[long_lag - 1] == 0) ? 1 : 0;
591 _M_p = 0;
592      }
593
594  template<typename _UIntType, size_t __w, size_t __s, size_t __r>
595    typename subtract_with_carry_engine<_UIntType, __w, __s, __r>::
596      result_type
597    subtract_with_carry_engine<_UIntType, __w, __s, __r>::
598    operator()()
599    {
600      // Derive short lag index from current index.
601      long __ps = _M_p - short_lag;
602      if (__ps < 0)
603 __ps += long_lag;
604
605      // Calculate new x(i) without overflow or division.
606      // NB: Thanks to the requirements for _UIntType, _M_x[_M_p] + _M_carry
607      // cannot overflow.
608      _UIntType __xi;
609      if (_M_x[__ps] >= _M_x[_M_p] + _M_carry)
610 {
611   __xi = _M_x[__ps] - _M_x[_M_p] - _M_carry;
612   _M_carry = 0;
613 }
614      else
615 {
616   __xi = (__detail::_Shift<_UIntType, __w>::__value
617   - _M_x[_M_p] - _M_carry + _M_x[__ps]);
618   _M_carry = 1;
619 }
620      _M_x[_M_p] = __xi;
621
622      // Adjust current index to loop around in ring buffer.
623      if (++_M_p >= long_lag)
624 _M_p = 0;
625
626      return __xi;
627    }
628
629  template<typename _UIntType, size_t __w, size_t __s, size_t __r,
630    typename _CharT, typename _Traits>
631    std::basic_ostream<_CharT, _Traits>&
632    operator<<(std::basic_ostream<_CharT, _Traits>& __os,
633        const subtract_with_carry_engine<_UIntType,
634 __w, __s, __r>& __x)
635    {
636      typedef std::basic_ostream<_CharT, _Traits>  __ostream_type;
637      typedef typename __ostream_type::ios_base    __ios_base;
638
639      const typename __ios_base::fmtflags __flags = __os.flags();
640      const _CharT __fill = __os.fill();
641      const _CharT __space = __os.widen(' ');
642      __os.flags(__ios_base::dec | __ios_base::fixed | __ios_base::left);
643      __os.fill(__space);
644
645      for (size_t __i = 0; __i < __r; ++__i)
646 __os << __x._M_x[__i] << __space;
647      __os << __x._M_carry << __space << __x._M_p;
648
649      __os.flags(__flags);
650      __os.fill(__fill);
651      return __os;
652    }
653
654  template<typename _UIntType, size_t __w, size_t __s, size_t __r,
655    typename _CharT, typename _Traits>
656    std::basic_istream<_CharT, _Traits>&
657    operator>>(std::basic_istream<_CharT, _Traits>& __is,
658        subtract_with_carry_engine<_UIntType, __w, __s, __r>& __x)
659    {
660      typedef std::basic_ostream<_CharT, _Traits>  __istream_type;
661      typedef typename __istream_type::ios_base    __ios_base;
662
663      const typename __ios_base::fmtflags __flags = __is.flags();
664      __is.flags(__ios_base::dec | __ios_base::skipws);
665
666      for (size_t __i = 0; __i < __r; ++__i)
667 __is >> __x._M_x[__i];
668      __is >> __x._M_carry;
669      __is >> __x._M_p;
670
671      __is.flags(__flags);
672      return __is;
673    }
674
675
676  template<typename _RandomNumberEngine, size_t __p, size_t __r>
677    constexpr size_t
678    discard_block_engine<_RandomNumberEngine, __p, __r>::block_size;
679
680  template<typename _RandomNumberEngine, size_t __p, size_t __r>
681    constexpr size_t
682    discard_block_engine<_RandomNumberEngine, __p, __r>::used_block;
683
684  template<typename _RandomNumberEngine, size_t __p, size_t __r>
685    typename discard_block_engine<_RandomNumberEngine,
686    __p, __r>::result_type
687    discard_block_engine<_RandomNumberEngine, __p, __r>::
688    operator()()
689    {
690      if (_M_n >= used_block)
691 {
692   _M_b.discard(block_size - _M_n);
693   _M_n = 0;
694 }
695      ++_M_n;
696      return _M_b();
697    }
698
699  template<typename _RandomNumberEngine, size_t __p, size_t __r,
700    typename _CharT, typename _Traits>
701    std::basic_ostream<_CharT, _Traits>&
702    operator<<(std::basic_ostream<_CharT, _Traits>& __os,
703        const discard_block_engine<_RandomNumberEngine,
704        __p, __r>& __x)
705    {
706      typedef std::basic_ostream<_CharT, _Traits>  __ostream_type;
707      typedef typename __ostream_type::ios_base    __ios_base;
708
709      const typename __ios_base::fmtflags __flags = __os.flags();
710      const _CharT __fill = __os.fill();
711      const _CharT __space = __os.widen(' ');
712      __os.flags(__ios_base::dec | __ios_base::fixed | __ios_base::left);
713      __os.fill(__space);
714
715      __os << __x.base() << __space << __x._M_n;
716
717      __os.flags(__flags);
718      __os.fill(__fill);
719      return __os;
720    }
721
722  template<typename _RandomNumberEngine, size_t __p, size_t __r,
723    typename _CharT, typename _Traits>
724    std::basic_istream<_CharT, _Traits>&
725    operator>>(std::basic_istream<_CharT, _Traits>& __is,
726        discard_block_engine<_RandomNumberEngine, __p, __r>& __x)
727    {
728      typedef std::basic_istream<_CharT, _Traits>  __istream_type;
729      typedef typename __istream_type::ios_base    __ios_base;
730
731      const typename __ios_base::fmtflags __flags = __is.flags();
732      __is.flags(__ios_base::dec | __ios_base::skipws);
733
734      __is >> __x._M_b >> __x._M_n;
735
736      __is.flags(__flags);
737      return __is;
738    }
739
740
741  template<typename _RandomNumberEngine, size_t __w, typename _UIntType>
742    typename independent_bits_engine<_RandomNumberEngine, __w, _UIntType>::
743      result_type
744    independent_bits_engine<_RandomNumberEngine, __w, _UIntType>::
745    operator()()
746    {
747      typedef typename _RandomNumberEngine::result_type _Eresult_type;
748      const _Eresult_type __r
749 = (_M_b.max() - _M_b.min() < std::numeric_limits<_Eresult_type>::max()
750    ? _M_b.max() - _M_b.min() + 1 : 0);
751      const unsigned __edig = std::numeric_limits<_Eresult_type>::digits;
752      const unsigned __m = __r ? std::__lg(__r) : __edig;
753
754      typedef typename std::common_type<_Eresult_type, result_type>::type
755 __ctype;
756      const unsigned __cdig = std::numeric_limits<__ctype>::digits;
757
758      unsigned __n, __n0;
759      __ctype __s0, __s1, __y0, __y1;
760
761      for (size_t __i = 0; __i < 2; ++__i)
762 {
763   __n = (__w + __m - 1) / __m + __i;
764   __n0 = __n - __w % __n;
765   const unsigned __w0 = __w / __n;  // __w0 <= __m
766
767   __s0 = 0;
768   __s1 = 0;
769   if (__w0 < __cdig)
770     {
771       __s0 = __ctype(1) << __w0;
772       __s1 = __s0 << 1;
773     }
774
775   __y0 = 0;
776   __y1 = 0;
777   if (__r)
778     {
779       __y0 = __s0 * (__r / __s0);
780       if (__s1)
781 __y1 = __s1 * (__r / __s1);
782
783       if (__r - __y0 <= __y0 / __n)
784 break;
785     }
786   else
787     break;
788 }
789
790      result_type __sum = 0;
791      for (size_t __k = 0; __k < __n0; ++__k)
792 {
793   __ctype __u;
794   do
795     __u = _M_b() - _M_b.min();
796   while (__y0 && __u >= __y0);
797   __sum = __s0 * __sum + (__s0 ? __u % __s0 : __u);
798 }
799      for (size_t __k = __n0; __k < __n; ++__k)
800 {
801   __ctype __u;
802   do
803     __u = _M_b() - _M_b.min();
804   while (__y1 && __u >= __y1);
805   __sum = __s1 * __sum + (__s1 ? __u % __s1 : __u);
806 }
807      return __sum;
808    }
809
810
811  template<typename _RandomNumberEngine, size_t __k>
812    constexpr size_t
813    shuffle_order_engine<_RandomNumberEngine, __k>::table_size;
814
815  template<typename _RandomNumberEngine, size_t __k>
816    typename shuffle_order_engine<_RandomNumberEngine, __k>::result_type
817    shuffle_order_engine<_RandomNumberEngine, __k>::
818    operator()()
819    {
820      size_t __j = __k * ((_M_y - _M_b.min())
821   / (_M_b.max() - _M_b.min() + 1.0L));
822      _M_y = _M_v[__j];
823      _M_v[__j] = _M_b();
824
825      return _M_y;
826    }
827
828  template<typename _RandomNumberEngine, size_t __k,
829    typename _CharT, typename _Traits>
830    std::basic_ostream<_CharT, _Traits>&
831    operator<<(std::basic_ostream<_CharT, _Traits>& __os,
832        const shuffle_order_engine<_RandomNumberEngine, __k>& __x)
833    {
834      typedef std::basic_ostream<_CharT, _Traits>  __ostream_type;
835      typedef typename __ostream_type::ios_base    __ios_base;
836
837      const typename __ios_base::fmtflags __flags = __os.flags();
838      const _CharT __fill = __os.fill();
839      const _CharT __space = __os.widen(' ');
840      __os.flags(__ios_base::dec | __ios_base::fixed | __ios_base::left);
841      __os.fill(__space);
842
843      __os << __x.base();
844      for (size_t __i = 0; __i < __k; ++__i)
845 __os << __space << __x._M_v[__i];
846      __os << __space << __x._M_y;
847
848      __os.flags(__flags);
849      __os.fill(__fill);
850      return __os;
851    }
852
853  template<typename _RandomNumberEngine, size_t __k,
854    typename _CharT, typename _Traits>
855    std::basic_istream<_CharT, _Traits>&
856    operator>>(std::basic_istream<_CharT, _Traits>& __is,
857        shuffle_order_engine<_RandomNumberEngine, __k>& __x)
858    {
859      typedef std::basic_istream<_CharT, _Traits>  __istream_type;
860      typedef typename __istream_type::ios_base    __ios_base;
861
862      const typename __ios_base::fmtflags __flags = __is.flags();
863      __is.flags(__ios_base::dec | __ios_base::skipws);
864
865      __is >> __x._M_b;
866      for (size_t __i = 0; __i < __k; ++__i)
867 __is >> __x._M_v[__i];
868      __is >> __x._M_y;
869
870      __is.flags(__flags);
871      return __is;
872    }
873
874
875  template<typename _IntType, typename _CharT, typename _Traits>
876    std::basic_ostream<_CharT, _Traits>&
877    operator<<(std::basic_ostream<_CharT, _Traits>& __os,
878        const uniform_int_distribution<_IntType>& __x)
879    {
880      typedef std::basic_ostream<_CharT, _Traits>  __ostream_type;
881      typedef typename __ostream_type::ios_base    __ios_base;
882
883      const typename __ios_base::fmtflags __flags = __os.flags();
884      const _CharT __fill = __os.fill();
885      const _CharT __space = __os.widen(' ');
886      __os.flags(__ios_base::scientific | __ios_base::left);
887      __os.fill(__space);
888
889      __os << __x.a() << __space << __x.b();
890
891      __os.flags(__flags);
892      __os.fill(__fill);
893      return __os;
894    }
895
896  template<typename _IntType, typename _CharT, typename _Traits>
897    std::basic_istream<_CharT, _Traits>&
898    operator>>(std::basic_istream<_CharT, _Traits>& __is,
899        uniform_int_distribution<_IntType>& __x)
900    {
901      typedef std::basic_istream<_CharT, _Traits>  __istream_type;
902      typedef typename __istream_type::ios_base    __ios_base;
903
904      const typename __ios_base::fmtflags __flags = __is.flags();
905      __is.flags(__ios_base::dec | __ios_base::skipws);
906
907      _IntType __a, __b;
908      __is >> __a >> __b;
909      __x.param(typename uniform_int_distribution<_IntType>::
910 param_type(__a, __b));
911
912      __is.flags(__flags);
913      return __is;
914    }
915
916
917  template<typename _RealType>
918    template<typename _ForwardIterator,
919      typename _UniformRandomNumberGenerator>
920      void
921      uniform_real_distribution<_RealType>::
922      __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
923       _UniformRandomNumberGenerator& __urng,
924       const param_type& __p)
925      {
926 __glibcxx_function_requires(_ForwardIteratorConcept<_ForwardIterator>)
927 __detail::_Adaptor<_UniformRandomNumberGenerator, result_type>
928   __aurng(__urng);
929 auto __range = __p.b() - __p.a();
930 while (__f != __t)
931   *__f++ = __aurng() * __range + __p.a();
932      }
933
934  template<typename _RealType, typename _CharT, typename _Traits>
935    std::basic_ostream<_CharT, _Traits>&
936    operator<<(std::basic_ostream<_CharT, _Traits>& __os,
937        const uniform_real_distribution<_RealType>& __x)
938    {
939      typedef std::basic_ostream<_CharT, _Traits>  __ostream_type;
940      typedef typename __ostream_type::ios_base    __ios_base;
941
942      const typename __ios_base::fmtflags __flags = __os.flags();
943      const _CharT __fill = __os.fill();
944      const std::streamsize __precision = __os.precision();
945      const _CharT __space = __os.widen(' ');
946      __os.flags(__ios_base::scientific | __ios_base::left);
947      __os.fill(__space);
948      __os.precision(std::numeric_limits<_RealType>::max_digits10);
949
950      __os << __x.a() << __space << __x.b();
951
952      __os.flags(__flags);
953      __os.fill(__fill);
954      __os.precision(__precision);
955      return __os;
956    }
957
958  template<typename _RealType, typename _CharT, typename _Traits>
959    std::basic_istream<_CharT, _Traits>&
960    operator>>(std::basic_istream<_CharT, _Traits>& __is,
961        uniform_real_distribution<_RealType>& __x)
962    {
963      typedef std::basic_istream<_CharT, _Traits>  __istream_type;
964      typedef typename __istream_type::ios_base    __ios_base;
965
966      const typename __ios_base::fmtflags __flags = __is.flags();
967      __is.flags(__ios_base::skipws);
968
969      _RealType __a, __b;
970      __is >> __a >> __b;
971      __x.param(typename uniform_real_distribution<_RealType>::
972 param_type(__a, __b));
973
974      __is.flags(__flags);
975      return __is;
976    }
977
978
979  template<typename _ForwardIterator,
980    typename _UniformRandomNumberGenerator>
981    void
982    std::bernoulli_distribution::
983    __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
984     _UniformRandomNumberGenerator& __urng,
985     const param_type& __p)
986    {
987      __glibcxx_function_requires(_ForwardIteratorConcept<_ForwardIterator>)
988      __detail::_Adaptor<_UniformRandomNumberGenerator, double>
989 __aurng(__urng);
990      auto __limit = __p.p() * (__aurng.max() - __aurng.min());
991
992      while (__f != __t)
993 *__f++ = (__aurng() - __aurng.min()) < __limit;
994    }
995
996  template<typename _CharT, typename _Traits>
997    std::basic_ostream<_CharT, _Traits>&
998    operator<<(std::basic_ostream<_CharT, _Traits>& __os,
999        const bernoulli_distribution& __x)
1000    {
1001      typedef std::basic_ostream<_CharT, _Traits>  __ostream_type;
1002      typedef typename __ostream_type::ios_base    __ios_base;
1003
1004      const typename __ios_base::fmtflags __flags = __os.flags();
1005      const _CharT __fill = __os.fill();
1006      const std::streamsize __precision = __os.precision();
1007      __os.flags(__ios_base::scientific | __ios_base::left);
1008      __os.fill(__os.widen(' '));
1009      __os.precision(std::numeric_limits<double>::max_digits10);
1010
1011      __os << __x.p();
1012
1013      __os.flags(__flags);
1014      __os.fill(__fill);
1015      __os.precision(__precision);
1016      return __os;
1017    }
1018
1019
1020  template<typename _IntType>
1021    template<typename _UniformRandomNumberGenerator>
1022      typename geometric_distribution<_IntType>::result_type
1023      geometric_distribution<_IntType>::
1024      operator()(_UniformRandomNumberGenerator& __urng,
1025  const param_type& __param)
1026      {
1027 // About the epsilon thing see this thread:
1028 // http://gcc.gnu.org/ml/gcc-patches/2006-10/msg00971.html
1029 const double __naf =
1030   (1 - std::numeric_limits<double>::epsilon()) / 2;
1031 // The largest _RealType convertible to _IntType.
1032 const double __thr =
1033   std::numeric_limits<_IntType>::max() + __naf;
1034 __detail::_Adaptor<_UniformRandomNumberGenerator, double>
1035   __aurng(__urng);
1036
1037 double __cand;
1038 do
1039   __cand = std::floor(std::log(1.0 - __aurng()) / __param._M_log_1_p);
1040 while (__cand >= __thr);
1041
1042 return result_type(__cand + __naf);
1043      }
1044
1045  template<typename _IntType>
1046    template<typename _ForwardIterator,
1047      typename _UniformRandomNumberGenerator>
1048      void
1049      geometric_distribution<_IntType>::
1050      __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
1051       _UniformRandomNumberGenerator& __urng,
1052       const param_type& __param)
1053      {
1054 __glibcxx_function_requires(_ForwardIteratorConcept<_ForwardIterator>)
1055 // About the epsilon thing see this thread:
1056 // http://gcc.gnu.org/ml/gcc-patches/2006-10/msg00971.html
1057 const double __naf =
1058   (1 - std::numeric_limits<double>::epsilon()) / 2;
1059 // The largest _RealType convertible to _IntType.
1060 const double __thr =
1061   std::numeric_limits<_IntType>::max() + __naf;
1062 __detail::_Adaptor<_UniformRandomNumberGenerator, double>
1063   __aurng(__urng);
1064
1065 while (__f != __t)
1066   {
1067     double __cand;
1068     do
1069       __cand = std::floor(std::log(1.0 - __aurng())
1070   / __param._M_log_1_p);
1071     while (__cand >= __thr);
1072
1073     *__f++ = __cand + __naf;
1074   }
1075      }
1076
1077  template<typename _IntType,
1078    typename _CharT, typename _Traits>
1079    std::basic_ostream<_CharT, _Traits>&
1080    operator<<(std::basic_ostream<_CharT, _Traits>& __os,
1081        const geometric_distribution<_IntType>& __x)
1082    {
1083      typedef std::basic_ostream<_CharT, _Traits>  __ostream_type;
1084      typedef typename __ostream_type::ios_base    __ios_base;
1085
1086      const typename __ios_base::fmtflags __flags = __os.flags();
1087      const _CharT __fill = __os.fill();
1088      const std::streamsize __precision = __os.precision();
1089      __os.flags(__ios_base::scientific | __ios_base::left);
1090      __os.fill(__os.widen(' '));
1091      __os.precision(std::numeric_limits<double>::max_digits10);
1092
1093      __os << __x.p();
1094
1095      __os.flags(__flags);
1096      __os.fill(__fill);
1097      __os.precision(__precision);
1098      return __os;
1099    }
1100
1101  template<typename _IntType,
1102    typename _CharT, typename _Traits>
1103    std::basic_istream<_CharT, _Traits>&
1104    operator>>(std::basic_istream<_CharT, _Traits>& __is,
1105        geometric_distribution<_IntType>& __x)
1106    {
1107      typedef std::basic_istream<_CharT, _Traits>  __istream_type;
1108      typedef typename __istream_type::ios_base    __ios_base;
1109
1110      const typename __ios_base::fmtflags __flags = __is.flags();
1111      __is.flags(__ios_base::skipws);
1112
1113      double __p;
1114      __is >> __p;
1115      __x.param(typename geometric_distribution<_IntType>::param_type(__p));
1116
1117      __is.flags(__flags);
1118      return __is;
1119    }
1120
1121  // This is Leger's algorithm, also in Devroye, Ch. X, Example 1.5.
1122  template<typename _IntType>
1123    template<typename _UniformRandomNumberGenerator>
1124      typename negative_binomial_distribution<_IntType>::result_type
1125      negative_binomial_distribution<_IntType>::
1126      operator()(_UniformRandomNumberGenerator& __urng)
1127      {
1128 const double __y = _M_gd(__urng);
1129
1130 // XXX Is the constructor too slow?
1131 std::poisson_distribution<result_type> __poisson(__y);
1132 return __poisson(__urng);
1133      }
1134
1135  template<typename _IntType>
1136    template<typename _UniformRandomNumberGenerator>
1137      typename negative_binomial_distribution<_IntType>::result_type
1138      negative_binomial_distribution<_IntType>::
1139      operator()(_UniformRandomNumberGenerator& __urng,
1140  const param_type& __p)
1141      {
1142 typedef typename std::gamma_distribution<double>::param_type
1143   param_type;
1144
1145 const double __y =
1146   _M_gd(__urng, param_type(__p.k(), (1.0 - __p.p()) / __p.p()));
1147
1148 std::poisson_distribution<result_type> __poisson(__y);
1149 return __poisson(__urng);
1150      }
1151
1152  template<typename _IntType>
1153    template<typename _ForwardIterator,
1154      typename _UniformRandomNumberGenerator>
1155      void
1156      negative_binomial_distribution<_IntType>::
1157      __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
1158       _UniformRandomNumberGenerator& __urng)
1159      {
1160 __glibcxx_function_requires(_ForwardIteratorConcept<_ForwardIterator>)
1161 while (__f != __t)
1162   {
1163     const double __y = _M_gd(__urng);
1164
1165     // XXX Is the constructor too slow?
1166     std::poisson_distribution<result_type> __poisson(__y);
1167     *__f++ = __poisson(__urng);
1168   }
1169      }
1170
1171  template<typename _IntType>
1172    template<typename _ForwardIterator,
1173      typename _UniformRandomNumberGenerator>
1174      void
1175      negative_binomial_distribution<_IntType>::
1176      __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
1177       _UniformRandomNumberGenerator& __urng,
1178       const param_type& __p)
1179      {
1180 __glibcxx_function_requires(_ForwardIteratorConcept<_ForwardIterator>)
1181 typename std::gamma_distribution<result_type>::param_type
1182   __p2(__p.k(), (1.0 - __p.p()) / __p.p());
1183
1184 while (__f != __t)
1185   {
1186     const double __y = _M_gd(__urng, __p2);
1187
1188     std::poisson_distribution<result_type> __poisson(__y);
1189     *__f++ = __poisson(__urng);
1190   }
1191      }
1192
1193  template<typename _IntType, typename _CharT, typename _Traits>
1194    std::basic_ostream<_CharT, _Traits>&
1195    operator<<(std::basic_ostream<_CharT, _Traits>& __os,
1196        const negative_binomial_distribution<_IntType>& __x)
1197    {
1198      typedef std::basic_ostream<_CharT, _Traits>  __ostream_type;
1199      typedef typename __ostream_type::ios_base    __ios_base;
1200
1201      const typename __ios_base::fmtflags __flags = __os.flags();
1202      const _CharT __fill = __os.fill();
1203      const std::streamsize __precision = __os.precision();
1204      const _CharT __space = __os.widen(' ');
1205      __os.flags(__ios_base::scientific | __ios_base::left);
1206      __os.fill(__os.widen(' '));
1207      __os.precision(std::numeric_limits<double>::max_digits10);
1208
1209      __os << __x.k() << __space << __x.p()
1210    << __space << __x._M_gd;
1211
1212      __os.flags(__flags);
1213      __os.fill(__fill);
1214      __os.precision(__precision);
1215      return __os;
1216    }
1217
1218  template<typename _IntType, typename _CharT, typename _Traits>
1219    std::basic_istream<_CharT, _Traits>&
1220    operator>>(std::basic_istream<_CharT, _Traits>& __is,
1221        negative_binomial_distribution<_IntType>& __x)
1222    {
1223      typedef std::basic_istream<_CharT, _Traits>  __istream_type;
1224      typedef typename __istream_type::ios_base    __ios_base;
1225
1226      const typename __ios_base::fmtflags __flags = __is.flags();
1227      __is.flags(__ios_base::skipws);
1228
1229      _IntType __k;
1230      double __p;
1231      __is >> __k >> __p >> __x._M_gd;
1232      __x.param(typename negative_binomial_distribution<_IntType>::
1233 param_type(__k, __p));
1234
1235      __is.flags(__flags);
1236      return __is;
1237    }
1238
1239
1240  template<typename _IntType>
1241    void
1242    poisson_distribution<_IntType>::param_type::
1243    _M_initialize()
1244    {
1245#if _GLIBCXX_USE_C99_MATH_TR1
1246      if (_M_mean >= 12)
1247 {
1248   const double __m = std::floor(_M_mean);
1249   _M_lm_thr = std::log(_M_mean);
1250   _M_lfm = std::lgamma(__m + 1);
1251   _M_sm = std::sqrt(__m);
1252
1253   const double __pi_4 = 0.7853981633974483096156608458198757L;
1254   const double __dx = std::sqrt(2 * __m * std::log(32 * __m
1255       / __pi_4));
1256   _M_d = std::round(std::max<double>(6.0, std::min(__m, __dx)));
1257   const double __cx = 2 * __m + _M_d;
1258   _M_scx = std::sqrt(__cx / 2);
1259   _M_1cx = 1 / __cx;
1260
1261   _M_c2b = std::sqrt(__pi_4 * __cx) * std::exp(_M_1cx);
1262   _M_cb = 2 * __cx * std::exp(-_M_d * _M_1cx * (1 + _M_d / 2))
1263 / _M_d;
1264 }
1265      else
1266#endif
1267 _M_lm_thr = std::exp(-_M_mean);
1268      }
1269
1270  /**
1271   * A rejection algorithm when mean >= 12 and a simple method based
1272   * upon the multiplication of uniform random variates otherwise.
1273   * NB: The former is available only if _GLIBCXX_USE_C99_MATH_TR1
1274   * is defined.
1275   *
1276   * Reference:
1277   * Devroye, L. Non-Uniform Random Variates Generation. Springer-Verlag,
1278   * New York, 1986, Ch. X, Sects. 3.3 & 3.4 (+ Errata!).
1279   */
1280  template<typename _IntType>
1281    template<typename _UniformRandomNumberGenerator>
1282      typename poisson_distribution<_IntType>::result_type
1283      poisson_distribution<_IntType>::
1284      operator()(_UniformRandomNumberGenerator& __urng,
1285  const param_type& __param)
1286      {
1287 __detail::_Adaptor<_UniformRandomNumberGenerator, double>
1288   __aurng(__urng);
1289#if _GLIBCXX_USE_C99_MATH_TR1
1290 if (__param.mean() >= 12)
1291   {
1292     double __x;
1293
1294     // See comments above...
1295     const double __naf =
1296       (1 - std::numeric_limits<double>::epsilon()) / 2;
1297     const double __thr =
1298       std::numeric_limits<_IntType>::max() + __naf;
1299
1300     const double __m = std::floor(__param.mean());
1301     // sqrt(pi / 2)
1302     const double __spi_2 = 1.2533141373155002512078826424055226L;
1303     const double __c1 = __param._M_sm * __spi_2;
1304     const double __c2 = __param._M_c2b + __c1;
1305     const double __c3 = __c2 + 1;
1306     const double __c4 = __c3 + 1;
1307     // e^(1 / 78)
1308     const double __e178 = 1.0129030479320018583185514777512983L;
1309     const double __c5 = __c4 + __e178;
1310     const double __c = __param._M_cb + __c5;
1311     const double __2cx = 2 * (2 * __m + __param._M_d);
1312
1313     bool __reject = true;
1314     do
1315       {
1316 const double __u = __c * __aurng();
1317 const double __e = -std::log(1.0 - __aurng());
1318
1319 double __w = 0.0;
1320
1321 if (__u <= __c1)
1322   {
1323     const double __n = _M_nd(__urng);
1324     const double __y = -std::abs(__n) * __param._M_sm - 1;
1325     __x = std::floor(__y);
1326     __w = -__n * __n / 2;
1327     if (__x < -__m)
1328       continue;
1329   }
1330 else if (__u <= __c2)
1331   {
1332     const double __n = _M_nd(__urng);
1333     const double __y = 1 + std::abs(__n) * __param._M_scx;
1334     __x = std::ceil(__y);
1335     __w = __y * (2 - __y) * __param._M_1cx;
1336     if (__x > __param._M_d)
1337       continue;
1338   }
1339 else if (__u <= __c3)
1340   // NB: This case not in the book, nor in the Errata,
1341   // but should be ok...
1342   __x = -1;
1343 else if (__u <= __c4)
1344   __x = 0;
1345 else if (__u <= __c5)
1346   __x = 1;
1347 else
1348   {
1349     const double __v = -std::log(1.0 - __aurng());
1350     const double __y = __param._M_d
1351      + __v * __2cx / __param._M_d;
1352     __x = std::ceil(__y);
1353     __w = -__param._M_d * __param._M_1cx * (1 + __y / 2);
1354   }
1355
1356 __reject = (__w - __e - __x * __param._M_lm_thr
1357     > __param._M_lfm - std::lgamma(__x + __m + 1));
1358
1359 __reject |= __x + __m >= __thr;
1360
1361       } while (__reject);
1362
1363     return result_type(__x + __m + __naf);
1364   }
1365 else
1366#endif
1367   {
1368     _IntType     __x = 0;
1369     double __prod = 1.0;
1370
1371     do
1372       {
1373 __prod *= __aurng();
1374 __x += 1;
1375       }
1376     while (__prod > __param._M_lm_thr);
1377
1378     return __x - 1;
1379   }
1380      }
1381
1382  template<typename _IntType>
1383    template<typename _ForwardIterator,
1384      typename _UniformRandomNumberGenerator>
1385      void
1386      poisson_distribution<_IntType>::
1387      __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
1388       _UniformRandomNumberGenerator& __urng,
1389       const param_type& __param)
1390      {
1391 __glibcxx_function_requires(_ForwardIteratorConcept<_ForwardIterator>)
1392 // We could duplicate everything from operator()...
1393 while (__f != __t)
1394   *__f++ = this->operator()(__urng, __param);
1395      }
1396
1397  template<typename _IntType,
1398    typename _CharT, typename _Traits>
1399    std::basic_ostream<_CharT, _Traits>&
1400    operator<<(std::basic_ostream<_CharT, _Traits>& __os,
1401        const poisson_distribution<_IntType>& __x)
1402    {
1403      typedef std::basic_ostream<_CharT, _Traits>  __ostream_type;
1404      typedef typename __ostream_type::ios_base    __ios_base;
1405
1406      const typename __ios_base::fmtflags __flags = __os.flags();
1407      const _CharT __fill = __os.fill();
1408      const std::streamsize __precision = __os.precision();
1409      const _CharT __space = __os.widen(' ');
1410      __os.flags(__ios_base::scientific | __ios_base::left);
1411      __os.fill(__space);
1412      __os.precision(std::numeric_limits<double>::max_digits10);
1413
1414      __os << __x.mean() << __space << __x._M_nd;
1415
1416      __os.flags(__flags);
1417      __os.fill(__fill);
1418      __os.precision(__precision);
1419      return __os;
1420    }
1421
1422  template<typename _IntType,
1423    typename _CharT, typename _Traits>
1424    std::basic_istream<_CharT, _Traits>&
1425    operator>>(std::basic_istream<_CharT, _Traits>& __is,
1426        poisson_distribution<_IntType>& __x)
1427    {
1428      typedef std::basic_istream<_CharT, _Traits>  __istream_type;
1429      typedef typename __istream_type::ios_base    __ios_base;
1430
1431      const typename __ios_base::fmtflags __flags = __is.flags();
1432      __is.flags(__ios_base::skipws);
1433
1434      double __mean;
1435      __is >> __mean >> __x._M_nd;
1436      __x.param(typename poisson_distribution<_IntType>::param_type(__mean));
1437
1438      __is.flags(__flags);
1439      return __is;
1440    }
1441
1442
1443  template<typename _IntType>
1444    void
1445    binomial_distribution<_IntType>::param_type::
1446    _M_initialize()
1447    {
1448      const double __p12 = _M_p <= 0.5 ? _M_p : 1.0 - _M_p;
1449
1450      _M_easy = true;
1451
1452#if _GLIBCXX_USE_C99_MATH_TR1
1453      if (_M_t * __p12 >= 8)
1454 {
1455   _M_easy = false;
1456   const double __np = std::floor(_M_t * __p12);
1457   const double __pa = __np / _M_t;
1458   const double __1p = 1 - __pa;
1459
1460   const double __pi_4 = 0.7853981633974483096156608458198757L;
1461   const double __d1x =
1462     std::sqrt(__np * __1p * std::log(32 * __np
1463      / (81 * __pi_4 * __1p)));
1464   _M_d1 = std::round(std::max<double>(1.0, __d1x));
1465   const double __d2x =
1466     std::sqrt(__np * __1p * std::log(32 * _M_t * __1p
1467      / (__pi_4 * __pa)));
1468   _M_d2 = std::round(std::max<double>(1.0, __d2x));
1469
1470   // sqrt(pi / 2)
1471   const double __spi_2 = 1.2533141373155002512078826424055226L;
1472   _M_s1 = std::sqrt(__np * __1p) * (1 + _M_d1 / (4 * __np));
1473   _M_s2 = std::sqrt(__np * __1p) * (1 + _M_d2 / (4 * _M_t * __1p));
1474   _M_c = 2 * _M_d1 / __np;
1475   _M_a1 = std::exp(_M_c) * _M_s1 * __spi_2;
1476   const double __a12 = _M_a1 + _M_s2 * __spi_2;
1477   const double __s1s = _M_s1 * _M_s1;
1478   _M_a123 = __a12 + (std::exp(_M_d1 / (_M_t * __1p))
1479      * 2 * __s1s / _M_d1
1480      * std::exp(-_M_d1 * _M_d1 / (2 * __s1s)));
1481   const double __s2s = _M_s2 * _M_s2;
1482   _M_s = (_M_a123 + 2 * __s2s / _M_d2
1483   * std::exp(-_M_d2 * _M_d2 / (2 * __s2s)));
1484   _M_lf = (std::lgamma(__np + 1)
1485    + std::lgamma(_M_t - __np + 1));
1486   _M_lp1p = std::log(__pa / __1p);
1487
1488   _M_q = -std::log(1 - (__p12 - __pa) / __1p);
1489 }
1490      else
1491#endif
1492 _M_q = -std::log(1 - __p12);
1493    }
1494
1495  template<typename _IntType>
1496    template<typename _UniformRandomNumberGenerator>
1497      typename binomial_distribution<_IntType>::result_type
1498      binomial_distribution<_IntType>::
1499      _M_waiting(_UniformRandomNumberGenerator& __urng,
1500  _IntType __t, double __q)
1501      {
1502 _IntType __x = 0;
1503 double __sum = 0.0;
1504 __detail::_Adaptor<_UniformRandomNumberGenerator, double>
1505   __aurng(__urng);
1506
1507 do
1508   {
1509     if (__t == __x)
1510       return __x;
1511     const double __e = -std::log(1.0 - __aurng());
1512     __sum += __e / (__t - __x);
1513     __x += 1;
1514   }
1515 while (__sum <= __q);
1516
1517 return __x - 1;
1518      }
1519
1520  /**
1521   * A rejection algorithm when t * p >= 8 and a simple waiting time
1522   * method - the second in the referenced book - otherwise.
1523   * NB: The former is available only if _GLIBCXX_USE_C99_MATH_TR1
1524   * is defined.
1525   *
1526   * Reference:
1527   * Devroye, L. Non-Uniform Random Variates Generation. Springer-Verlag,
1528   * New York, 1986, Ch. X, Sect. 4 (+ Errata!).
1529   */
1530  template<typename _IntType>
1531    template<typename _UniformRandomNumberGenerator>
1532      typename binomial_distribution<_IntType>::result_type
1533      binomial_distribution<_IntType>::
1534      operator()(_UniformRandomNumberGenerator& __urng,
1535  const param_type& __param)
1536      {
1537 result_type __ret;
1538 const _IntType __t = __param.t();
1539 const double __p = __param.p();
1540 const double __p12 = __p <= 0.5 ? __p : 1.0 - __p;
1541 __detail::_Adaptor<_UniformRandomNumberGenerator, double>
1542   __aurng(__urng);
1543
1544#if _GLIBCXX_USE_C99_MATH_TR1
1545 if (!__param._M_easy)
1546   {
1547     double __x;
1548
1549     // See comments above...
1550     const double __naf =
1551       (1 - std::numeric_limits<double>::epsilon()) / 2;
1552     const double __thr =
1553       std::numeric_limits<_IntType>::max() + __naf;
1554
1555     const double __np = std::floor(__t * __p12);
1556
1557     // sqrt(pi / 2)
1558     const double __spi_2 = 1.2533141373155002512078826424055226L;
1559     const double __a1 = __param._M_a1;
1560     const double __a12 = __a1 + __param._M_s2 * __spi_2;
1561     const double __a123 = __param._M_a123;
1562     const double __s1s = __param._M_s1 * __param._M_s1;
1563     const double __s2s = __param._M_s2 * __param._M_s2;
1564
1565     bool __reject;
1566     do
1567       {
1568 const double __u = __param._M_s * __aurng();
1569
1570 double __v;
1571
1572 if (__u <= __a1)
1573   {
1574     const double __n = _M_nd(__urng);
1575     const double __y = __param._M_s1 * std::abs(__n);
1576     __reject = __y >= __param._M_d1;
1577     if (!__reject)
1578       {
1579 const double __e = -std::log(1.0 - __aurng());
1580 __x = std::floor(__y);
1581 __v = -__e - __n * __n / 2 + __param._M_c;
1582       }
1583   }
1584 else if (__u <= __a12)
1585   {
1586     const double __n = _M_nd(__urng);
1587     const double __y = __param._M_s2 * std::abs(__n);
1588     __reject = __y >= __param._M_d2;
1589     if (!__reject)
1590       {
1591 const double __e = -std::log(1.0 - __aurng());
1592 __x = std::floor(-__y);
1593 __v = -__e - __n * __n / 2;
1594       }
1595   }
1596 else if (__u <= __a123)
1597   {
1598     const double __e1 = -std::log(1.0 - __aurng());
1599     const double __e2 = -std::log(1.0 - __aurng());
1600
1601     const double __y = __param._M_d1
1602      + 2 * __s1s * __e1 / __param._M_d1;
1603     __x = std::floor(__y);
1604     __v = (-__e2 + __param._M_d1 * (1 / (__t - __np)
1605     -__y / (2 * __s1s)));
1606     __reject = false;
1607   }
1608 else
1609   {
1610     const double __e1 = -std::log(1.0 - __aurng());
1611     const double __e2 = -std::log(1.0 - __aurng());
1612
1613     const double __y = __param._M_d2
1614      + 2 * __s2s * __e1 / __param._M_d2;
1615     __x = std::floor(-__y);
1616     __v = -__e2 - __param._M_d2 * __y / (2 * __s2s);
1617     __reject = false;
1618   }
1619
1620 __reject = __reject || __x < -__np || __x > __t - __np;
1621 if (!__reject)
1622   {
1623     const double __lfx =
1624       std::lgamma(__np + __x + 1)
1625       + std::lgamma(__t - (__np + __x) + 1);
1626     __reject = __v > __param._M_lf - __lfx
1627      + __x * __param._M_lp1p;
1628   }
1629
1630 __reject |= __x + __np >= __thr;
1631       }
1632     while (__reject);
1633
1634     __x += __np + __naf;
1635
1636     const _IntType __z = _M_waiting(__urng, __t - _IntType(__x),
1637     __param._M_q);
1638     __ret = _IntType(__x) + __z;
1639   }
1640 else
1641#endif
1642   __ret = _M_waiting(__urng, __t, __param._M_q);
1643
1644 if (__p12 != __p)
1645   __ret = __t - __ret;
1646 return __ret;
1647      }
1648
1649  template<typename _IntType>
1650    template<typename _ForwardIterator,
1651      typename _UniformRandomNumberGenerator>
1652      void
1653      binomial_distribution<_IntType>::
1654      __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
1655       _UniformRandomNumberGenerator& __urng,
1656       const param_type& __param)
1657      {
1658 __glibcxx_function_requires(_ForwardIteratorConcept<_ForwardIterator>)
1659 // We could duplicate everything from operator()...
1660 while (__f != __t)
1661   *__f++ = this->operator()(__urng, __param);
1662      }
1663
1664  template<typename _IntType,
1665    typename _CharT, typename _Traits>
1666    std::basic_ostream<_CharT, _Traits>&
1667    operator<<(std::basic_ostream<_CharT, _Traits>& __os,
1668        const binomial_distribution<_IntType>& __x)
1669    {
1670      typedef std::basic_ostream<_CharT, _Traits>  __ostream_type;
1671      typedef typename __ostream_type::ios_base    __ios_base;
1672
1673      const typename __ios_base::fmtflags __flags = __os.flags();
1674      const _CharT __fill = __os.fill();
1675      const std::streamsize __precision = __os.precision();
1676      const _CharT __space = __os.widen(' ');
1677      __os.flags(__ios_base::scientific | __ios_base::left);
1678      __os.fill(__space);
1679      __os.precision(std::numeric_limits<double>::max_digits10);
1680
1681      __os << __x.t() << __space << __x.p()
1682    << __space << __x._M_nd;
1683
1684      __os.flags(__flags);
1685      __os.fill(__fill);
1686      __os.precision(__precision);
1687      return __os;
1688    }
1689
1690  template<typename _IntType,
1691    typename _CharT, typename _Traits>
1692    std::basic_istream<_CharT, _Traits>&
1693    operator>>(std::basic_istream<_CharT, _Traits>& __is,
1694        binomial_distribution<_IntType>& __x)
1695    {
1696      typedef std::basic_istream<_CharT, _Traits>  __istream_type;
1697      typedef typename __istream_type::ios_base    __ios_base;
1698
1699      const typename __ios_base::fmtflags __flags = __is.flags();
1700      __is.flags(__ios_base::dec | __ios_base::skipws);
1701
1702      _IntType __t;
1703      double __p;
1704      __is >> __t >> __p >> __x._M_nd;
1705      __x.param(typename binomial_distribution<_IntType>::
1706 param_type(__t, __p));
1707
1708      __is.flags(__flags);
1709      return __is;
1710    }
1711
1712
1713  template<typename _RealType>
1714    template<typename _ForwardIterator,
1715      typename _UniformRandomNumberGenerator>
1716      void
1717      std::exponential_distribution<_RealType>::
1718      __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
1719       _UniformRandomNumberGenerator& __urng,
1720       const param_type& __p)
1721      {
1722 __glibcxx_function_requires(_ForwardIteratorConcept<_ForwardIterator>)
1723 __detail::_Adaptor<_UniformRandomNumberGenerator, result_type>
1724   __aurng(__urng);
1725 while (__f != __t)
1726   *__f++ = -std::log(result_type(1) - __aurng()) / __p.lambda();
1727      }
1728
1729  template<typename _RealType, typename _CharT, typename _Traits>
1730    std::basic_ostream<_CharT, _Traits>&
1731    operator<<(std::basic_ostream<_CharT, _Traits>& __os,
1732        const exponential_distribution<_RealType>& __x)
1733    {
1734      typedef std::basic_ostream<_CharT, _Traits>  __ostream_type;
1735      typedef typename __ostream_type::ios_base    __ios_base;
1736
1737      const typename __ios_base::fmtflags __flags = __os.flags();
1738      const _CharT __fill = __os.fill();
1739      const std::streamsize __precision = __os.precision();
1740      __os.flags(__ios_base::scientific | __ios_base::left);
1741      __os.fill(__os.widen(' '));
1742      __os.precision(std::numeric_limits<_RealType>::max_digits10);
1743
1744      __os << __x.lambda();
1745
1746      __os.flags(__flags);
1747      __os.fill(__fill);
1748      __os.precision(__precision);
1749      return __os;
1750    }
1751
1752  template<typename _RealType, typename _CharT, typename _Traits>
1753    std::basic_istream<_CharT, _Traits>&
1754    operator>>(std::basic_istream<_CharT, _Traits>& __is,
1755        exponential_distribution<_RealType>& __x)
1756    {
1757      typedef std::basic_istream<_CharT, _Traits>  __istream_type;
1758      typedef typename __istream_type::ios_base    __ios_base;
1759
1760      const typename __ios_base::fmtflags __flags = __is.flags();
1761      __is.flags(__ios_base::dec | __ios_base::skipws);
1762
1763      _RealType __lambda;
1764      __is >> __lambda;
1765      __x.param(typename exponential_distribution<_RealType>::
1766 param_type(__lambda));
1767
1768      __is.flags(__flags);
1769      return __is;
1770    }
1771
1772
1773  /**
1774   * Polar method due to Marsaglia.
1775   *
1776   * Devroye, L. Non-Uniform Random Variates Generation. Springer-Verlag,
1777   * New York, 1986, Ch. V, Sect. 4.4.
1778   */
1779  template<typename _RealType>
1780    template<typename _UniformRandomNumberGenerator>
1781      typename normal_distribution<_RealType>::result_type
1782      normal_distribution<_RealType>::
1783      operator()(_UniformRandomNumberGenerator& __urng,
1784  const param_type& __param)
1785      {
1786 result_type __ret;
1787 __detail::_Adaptor<_UniformRandomNumberGenerator, result_type>
1788   __aurng(__urng);
1789
1790 if (_M_saved_available)
1791   {
1792     _M_saved_available = false;
1793     __ret = _M_saved;
1794   }
1795 else
1796   {
1797     result_type __x, __y, __r2;
1798     do
1799       {
1800 __x = result_type(2.0) * __aurng() - 1.0;
1801 __y = result_type(2.0) * __aurng() - 1.0;
1802 __r2 = __x * __x + __y * __y;
1803       }
1804     while (__r2 > 1.0 || __r2 == 0.0);
1805
1806     const result_type __mult = std::sqrt(-2 * std::log(__r2) / __r2);
1807     _M_saved = __x * __mult;
1808     _M_saved_available = true;
1809     __ret = __y * __mult;
1810   }
1811
1812 __ret = __ret * __param.stddev() + __param.mean();
1813 return __ret;
1814      }
1815
1816  template<typename _RealType>
1817    template<typename _ForwardIterator,
1818      typename _UniformRandomNumberGenerator>
1819      void
1820      normal_distribution<_RealType>::
1821      __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
1822       _UniformRandomNumberGenerator& __urng,
1823       const param_type& __param)
1824      {
1825 __glibcxx_function_requires(_ForwardIteratorConcept<_ForwardIterator>)
1826
1827 if (__f == __t)
1828   return;
1829
1830 if (_M_saved_available)
1831   {
1832     _M_saved_available = false;
1833     *__f++ = _M_saved * __param.stddev() + __param.mean();
1834
1835     if (__f == __t)
1836       return;
1837   }
1838
1839 __detail::_Adaptor<_UniformRandomNumberGenerator, result_type>
1840   __aurng(__urng);
1841
1842 while (__f + 1 < __t)
1843   {
1844     result_type __x, __y, __r2;
1845     do
1846       {
1847 __x = result_type(2.0) * __aurng() - 1.0;
1848 __y = result_type(2.0) * __aurng() - 1.0;
1849 __r2 = __x * __x + __y * __y;
1850       }
1851     while (__r2 > 1.0 || __r2 == 0.0);
1852
1853     const result_type __mult = std::sqrt(-2 * std::log(__r2) / __r2);
1854     *__f++ = __y * __mult * __param.stddev() + __param.mean();
1855     *__f++ = __x * __mult * __param.stddev() + __param.mean();
1856   }
1857
1858 if (__f != __t)
1859   {
1860     result_type __x, __y, __r2;
1861     do
1862       {
1863 __x = result_type(2.0) * __aurng() - 1.0;
1864 __y = result_type(2.0) * __aurng() - 1.0;
1865 __r2 = __x * __x + __y * __y;
1866       }
1867     while (__r2 > 1.0 || __r2 == 0.0);
1868
1869     const result_type __mult = std::sqrt(-2 * std::log(__r2) / __r2);
1870     _M_saved = __x * __mult;
1871     _M_saved_available = true;
1872     *__f = __y * __mult * __param.stddev() + __param.mean();
1873   }
1874      }
1875
1876  template<typename _RealType>
1877    bool
1878    operator==(const std::normal_distribution<_RealType>& __d1,
1879        const std::normal_distribution<_RealType>& __d2)
1880    {
1881      if (__d1._M_param == __d2._M_param
1882   && __d1._M_saved_available == __d2._M_saved_available)
1883 {
1884   if (__d1._M_saved_available
1885       && __d1._M_saved == __d2._M_saved)
1886     return true;
1887   else if(!__d1._M_saved_available)
1888     return true;
1889   else
1890     return false;
1891 }
1892      else
1893 return false;
1894    }
1895
1896  template<typename _RealType, typename _CharT, typename _Traits>
1897    std::basic_ostream<_CharT, _Traits>&
1898    operator<<(std::basic_ostream<_CharT, _Traits>& __os,
1899        const normal_distribution<_RealType>& __x)
1900    {
1901      typedef std::basic_ostream<_CharT, _Traits>  __ostream_type;
1902      typedef typename __ostream_type::ios_base    __ios_base;
1903
1904      const typename __ios_base::fmtflags __flags = __os.flags();
1905      const _CharT __fill = __os.fill();
1906      const std::streamsize __precision = __os.precision();
1907      const _CharT __space = __os.widen(' ');
1908      __os.flags(__ios_base::scientific | __ios_base::left);
1909      __os.fill(__space);
1910      __os.precision(std::numeric_limits<_RealType>::max_digits10);
1911
1912      __os << __x.mean() << __space << __x.stddev()
1913    << __space << __x._M_saved_available;
1914      if (__x._M_saved_available)
1915 __os << __space << __x._M_saved;
1916
1917      __os.flags(__flags);
1918      __os.fill(__fill);
1919      __os.precision(__precision);
1920      return __os;
1921    }
1922
1923  template<typename _RealType, typename _CharT, typename _Traits>
1924    std::basic_istream<_CharT, _Traits>&
1925    operator>>(std::basic_istream<_CharT, _Traits>& __is,
1926        normal_distribution<_RealType>& __x)
1927    {
1928      typedef std::basic_istream<_CharT, _Traits>  __istream_type;
1929      typedef typename __istream_type::ios_base    __ios_base;
1930
1931      const typename __ios_base::fmtflags __flags = __is.flags();
1932      __is.flags(__ios_base::dec | __ios_base::skipws);
1933
1934      double __mean, __stddev;
1935      __is >> __mean >> __stddev
1936    >> __x._M_saved_available;
1937      if (__x._M_saved_available)
1938 __is >> __x._M_saved;
1939      __x.param(typename normal_distribution<_RealType>::
1940 param_type(__mean, __stddev));
1941
1942      __is.flags(__flags);
1943      return __is;
1944    }
1945
1946
1947  template<typename _RealType>
1948    template<typename _ForwardIterator,
1949      typename _UniformRandomNumberGenerator>
1950      void
1951      lognormal_distribution<_RealType>::
1952      __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
1953       _UniformRandomNumberGenerator& __urng,
1954       const param_type& __p)
1955      {
1956 __glibcxx_function_requires(_ForwardIteratorConcept<_ForwardIterator>)
1957   while (__f != __t)
1958     *__f++ = std::exp(__p.s() * _M_nd(__urng) + __p.m());
1959      }
1960
1961  template<typename _RealType, typename _CharT, typename _Traits>
1962    std::basic_ostream<_CharT, _Traits>&
1963    operator<<(std::basic_ostream<_CharT, _Traits>& __os,
1964        const lognormal_distribution<_RealType>& __x)
1965    {
1966      typedef std::basic_ostream<_CharT, _Traits>  __ostream_type;
1967      typedef typename __ostream_type::ios_base    __ios_base;
1968
1969      const typename __ios_base::fmtflags __flags = __os.flags();
1970      const _CharT __fill = __os.fill();
1971      const std::streamsize __precision = __os.precision();
1972      const _CharT __space = __os.widen(' ');
1973      __os.flags(__ios_base::scientific | __ios_base::left);
1974      __os.fill(__space);
1975      __os.precision(std::numeric_limits<_RealType>::max_digits10);
1976
1977      __os << __x.m() << __space << __x.s()
1978    << __space << __x._M_nd;
1979
1980      __os.flags(__flags);
1981      __os.fill(__fill);
1982      __os.precision(__precision);
1983      return __os;
1984    }
1985
1986  template<typename _RealType, typename _CharT, typename _Traits>
1987    std::basic_istream<_CharT, _Traits>&
1988    operator>>(std::basic_istream<_CharT, _Traits>& __is,
1989        lognormal_distribution<_RealType>& __x)
1990    {
1991      typedef std::basic_istream<_CharT, _Traits>  __istream_type;
1992      typedef typename __istream_type::ios_base    __ios_base;
1993
1994      const typename __ios_base::fmtflags __flags = __is.flags();
1995      __is.flags(__ios_base::dec | __ios_base::skipws);
1996
1997      _RealType __m, __s;
1998      __is >> __m >> __s >> __x._M_nd;
1999      __x.param(typename lognormal_distribution<_RealType>::
2000 param_type(__m, __s));
2001
2002      __is.flags(__flags);
2003      return __is;
2004    }
2005
2006  template<typename _RealType>
2007    template<typename _ForwardIterator,
2008      typename _UniformRandomNumberGenerator>
2009      void
2010      std::chi_squared_distribution<_RealType>::
2011      __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
2012       _UniformRandomNumberGenerator& __urng)
2013      {
2014 __glibcxx_function_requires(_ForwardIteratorConcept<_ForwardIterator>)
2015 while (__f != __t)
2016   *__f++ = 2 * _M_gd(__urng);
2017      }
2018
2019  template<typename _RealType>
2020    template<typename _ForwardIterator,
2021      typename _UniformRandomNumberGenerator>
2022      void
2023      std::chi_squared_distribution<_RealType>::
2024      __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
2025       _UniformRandomNumberGenerator& __urng,
2026       const typename
2027       std::gamma_distribution<result_type>::param_type& __p)
2028      {
2029 __glibcxx_function_requires(_ForwardIteratorConcept<_ForwardIterator>)
2030 while (__f != __t)
2031   *__f++ = 2 * _M_gd(__urng, __p);
2032      }
2033
2034  template<typename _RealType, typename _CharT, typename _Traits>
2035    std::basic_ostream<_CharT, _Traits>&
2036    operator<<(std::basic_ostream<_CharT, _Traits>& __os,
2037        const chi_squared_distribution<_RealType>& __x)
2038    {
2039      typedef std::basic_ostream<_CharT, _Traits>  __ostream_type;
2040      typedef typename __ostream_type::ios_base    __ios_base;
2041
2042      const typename __ios_base::fmtflags __flags = __os.flags();
2043      const _CharT __fill = __os.fill();
2044      const std::streamsize __precision = __os.precision();
2045      const _CharT __space = __os.widen(' ');
2046      __os.flags(__ios_base::scientific | __ios_base::left);
2047      __os.fill(__space);
2048      __os.precision(std::numeric_limits<_RealType>::max_digits10);
2049
2050      __os << __x.n() << __space << __x._M_gd;
2051
2052      __os.flags(__flags);
2053      __os.fill(__fill);
2054      __os.precision(__precision);
2055      return __os;
2056    }
2057
2058  template<typename _RealType, typename _CharT, typename _Traits>
2059    std::basic_istream<_CharT, _Traits>&
2060    operator>>(std::basic_istream<_CharT, _Traits>& __is,
2061        chi_squared_distribution<_RealType>& __x)
2062    {
2063      typedef std::basic_istream<_CharT, _Traits>  __istream_type;
2064      typedef typename __istream_type::ios_base    __ios_base;
2065
2066      const typename __ios_base::fmtflags __flags = __is.flags();
2067      __is.flags(__ios_base::dec | __ios_base::skipws);
2068
2069      _RealType __n;
2070      __is >> __n >> __x._M_gd;
2071      __x.param(typename chi_squared_distribution<_RealType>::
2072 param_type(__n));
2073
2074      __is.flags(__flags);
2075      return __is;
2076    }
2077
2078
2079  template<typename _RealType>
2080    template<typename _UniformRandomNumberGenerator>
2081      typename cauchy_distribution<_RealType>::result_type
2082      cauchy_distribution<_RealType>::
2083      operator()(_UniformRandomNumberGenerator& __urng,
2084  const param_type& __p)
2085      {
2086 __detail::_Adaptor<_UniformRandomNumberGenerator, result_type>
2087   __aurng(__urng);
2088 _RealType __u;
2089 do
2090   __u = __aurng();
2091 while (__u == 0.5);
2092
2093 const _RealType __pi = 3.1415926535897932384626433832795029L;
2094 return __p.a() + __p.b() * std::tan(__pi * __u);
2095      }
2096
2097  template<typename _RealType>
2098    template<typename _ForwardIterator,
2099      typename _UniformRandomNumberGenerator>
2100      void
2101      cauchy_distribution<_RealType>::
2102      __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
2103       _UniformRandomNumberGenerator& __urng,
2104       const param_type& __p)
2105      {
2106 __glibcxx_function_requires(_ForwardIteratorConcept<_ForwardIterator>)
2107 const _RealType __pi = 3.1415926535897932384626433832795029L;
2108 __detail::_Adaptor<_UniformRandomNumberGenerator, result_type>
2109   __aurng(__urng);
2110 while (__f != __t)
2111   {
2112     _RealType __u;
2113     do
2114       __u = __aurng();
2115     while (__u == 0.5);
2116
2117     *__f++ = __p.a() + __p.b() * std::tan(__pi * __u);
2118   }
2119      }
2120
2121  template<typename _RealType, typename _CharT, typename _Traits>
2122    std::basic_ostream<_CharT, _Traits>&
2123    operator<<(std::basic_ostream<_CharT, _Traits>& __os,
2124        const cauchy_distribution<_RealType>& __x)
2125    {
2126      typedef std::basic_ostream<_CharT, _Traits>  __ostream_type;
2127      typedef typename __ostream_type::ios_base    __ios_base;
2128
2129      const typename __ios_base::fmtflags __flags = __os.flags();
2130      const _CharT __fill = __os.fill();
2131      const std::streamsize __precision = __os.precision();
2132      const _CharT __space = __os.widen(' ');
2133      __os.flags(__ios_base::scientific | __ios_base::left);
2134      __os.fill(__space);
2135      __os.precision(std::numeric_limits<_RealType>::max_digits10);
2136
2137      __os << __x.a() << __space << __x.b();
2138
2139      __os.flags(__flags);
2140      __os.fill(__fill);
2141      __os.precision(__precision);
2142      return __os;
2143    }
2144
2145  template<typename _RealType, typename _CharT, typename _Traits>
2146    std::basic_istream<_CharT, _Traits>&
2147    operator>>(std::basic_istream<_CharT, _Traits>& __is,
2148        cauchy_distribution<_RealType>& __x)
2149    {
2150      typedef std::basic_istream<_CharT, _Traits>  __istream_type;
2151      typedef typename __istream_type::ios_base    __ios_base;
2152
2153      const typename __ios_base::fmtflags __flags = __is.flags();
2154      __is.flags(__ios_base::dec | __ios_base::skipws);
2155
2156      _RealType __a, __b;
2157      __is >> __a >> __b;
2158      __x.param(typename cauchy_distribution<_RealType>::
2159 param_type(__a, __b));
2160
2161      __is.flags(__flags);
2162      return __is;
2163    }
2164
2165
2166  template<typename _RealType>
2167    template<typename _ForwardIterator,
2168      typename _UniformRandomNumberGenerator>
2169      void
2170      std::fisher_f_distribution<_RealType>::
2171      __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
2172       _UniformRandomNumberGenerator& __urng)
2173      {
2174 __glibcxx_function_requires(_ForwardIteratorConcept<_ForwardIterator>)
2175 while (__f != __t)
2176   *__f++ = ((_M_gd_x(__urng) * n()) / (_M_gd_y(__urng) * m()));
2177      }
2178
2179  template<typename _RealType>
2180    template<typename _ForwardIterator,
2181      typename _UniformRandomNumberGenerator>
2182      void
2183      std::fisher_f_distribution<_RealType>::
2184      __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
2185       _UniformRandomNumberGenerator& __urng,
2186       const param_type& __p)
2187      {
2188 __glibcxx_function_requires(_ForwardIteratorConcept<_ForwardIterator>)
2189 typedef typename std::gamma_distribution<result_type>::param_type
2190   param_type;
2191 param_type __p1(__p.m() / 2);
2192 param_type __p2(__p.n() / 2);
2193 while (__f != __t)
2194   *__f++ = ((_M_gd_x(__urng, __p1) * n())
2195     / (_M_gd_y(__urng, __p2) * m()));
2196      }
2197
2198  template<typename _RealType, typename _CharT, typename _Traits>
2199    std::basic_ostream<_CharT, _Traits>&
2200    operator<<(std::basic_ostream<_CharT, _Traits>& __os,
2201        const fisher_f_distribution<_RealType>& __x)
2202    {
2203      typedef std::basic_ostream<_CharT, _Traits>  __ostream_type;
2204      typedef typename __ostream_type::ios_base    __ios_base;
2205
2206      const typename __ios_base::fmtflags __flags = __os.flags();
2207      const _CharT __fill = __os.fill();
2208      const std::streamsize __precision = __os.precision();
2209      const _CharT __space = __os.widen(' ');
2210      __os.flags(__ios_base::scientific | __ios_base::left);
2211      __os.fill(__space);
2212      __os.precision(std::numeric_limits<_RealType>::max_digits10);
2213
2214      __os << __x.m() << __space << __x.n()
2215    << __space << __x._M_gd_x << __space << __x._M_gd_y;
2216
2217      __os.flags(__flags);
2218      __os.fill(__fill);
2219      __os.precision(__precision);
2220      return __os;
2221    }
2222
2223  template<typename _RealType, typename _CharT, typename _Traits>
2224    std::basic_istream<_CharT, _Traits>&
2225    operator>>(std::basic_istream<_CharT, _Traits>& __is,
2226        fisher_f_distribution<_RealType>& __x)
2227    {
2228      typedef std::basic_istream<_CharT, _Traits>  __istream_type;
2229      typedef typename __istream_type::ios_base    __ios_base;
2230
2231      const typename __ios_base::fmtflags __flags = __is.flags();
2232      __is.flags(__ios_base::dec | __ios_base::skipws);
2233
2234      _RealType __m, __n;
2235      __is >> __m >> __n >> __x._M_gd_x >> __x._M_gd_y;
2236      __x.param(typename fisher_f_distribution<_RealType>::
2237 param_type(__m, __n));
2238
2239      __is.flags(__flags);
2240      return __is;
2241    }
2242
2243
2244  template<typename _RealType>
2245    template<typename _ForwardIterator,
2246      typename _UniformRandomNumberGenerator>
2247      void
2248      std::student_t_distribution<_RealType>::
2249      __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
2250       _UniformRandomNumberGenerator& __urng)
2251      {
2252 __glibcxx_function_requires(_ForwardIteratorConcept<_ForwardIterator>)
2253 while (__f != __t)
2254   *__f++ = _M_nd(__urng) * std::sqrt(n() / _M_gd(__urng));
2255      }
2256
2257  template<typename _RealType>
2258    template<typename _ForwardIterator,
2259      typename _UniformRandomNumberGenerator>
2260      void
2261      std::student_t_distribution<_RealType>::
2262      __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
2263       _UniformRandomNumberGenerator& __urng,
2264       const param_type& __p)
2265      {
2266 __glibcxx_function_requires(_ForwardIteratorConcept<_ForwardIterator>)
2267 typename std::gamma_distribution<result_type>::param_type
2268   __p2(__p.n() / 2, 2);
2269 while (__f != __t)
2270   *__f++ =  _M_nd(__urng) * std::sqrt(__p.n() / _M_gd(__urng, __p2));
2271      }
2272
2273  template<typename _RealType, typename _CharT, typename _Traits>
2274    std::basic_ostream<_CharT, _Traits>&
2275    operator<<(std::basic_ostream<_CharT, _Traits>& __os,
2276        const student_t_distribution<_RealType>& __x)
2277    {
2278      typedef std::basic_ostream<_CharT, _Traits>  __ostream_type;
2279      typedef typename __ostream_type::ios_base    __ios_base;
2280
2281      const typename __ios_base::fmtflags __flags = __os.flags();
2282      const _CharT __fill = __os.fill();
2283      const std::streamsize __precision = __os.precision();
2284      const _CharT __space = __os.widen(' ');
2285      __os.flags(__ios_base::scientific | __ios_base::left);
2286      __os.fill(__space);
2287      __os.precision(std::numeric_limits<_RealType>::max_digits10);
2288
2289      __os << __x.n() << __space << __x._M_nd << __space << __x._M_gd;
2290
2291      __os.flags(__flags);
2292      __os.fill(__fill);
2293      __os.precision(__precision);
2294      return __os;
2295    }
2296
2297  template<typename _RealType, typename _CharT, typename _Traits>
2298    std::basic_istream<_CharT, _Traits>&
2299    operator>>(std::basic_istream<_CharT, _Traits>& __is,
2300        student_t_distribution<_RealType>& __x)
2301    {
2302      typedef std::basic_istream<_CharT, _Traits>  __istream_type;
2303      typedef typename __istream_type::ios_base    __ios_base;
2304
2305      const typename __ios_base::fmtflags __flags = __is.flags();
2306      __is.flags(__ios_base::dec | __ios_base::skipws);
2307
2308      _RealType __n;
2309      __is >> __n >> __x._M_nd >> __x._M_gd;
2310      __x.param(typename student_t_distribution<_RealType>::param_type(__n));
2311
2312      __is.flags(__flags);
2313      return __is;
2314    }
2315
2316
2317  template<typename _RealType>
2318    void
2319    gamma_distribution<_RealType>::param_type::
2320    _M_initialize()
2321    {
2322      _M_malpha = _M_alpha < 1.0 ? _M_alpha + _RealType(1.0) : _M_alpha;
2323
2324      const _RealType __a1 = _M_malpha - _RealType(1.0) / _RealType(3.0);
2325      _M_a2 = _RealType(1.0) / std::sqrt(_RealType(9.0) * __a1);
2326    }
2327
2328  /**
2329   * Marsaglia, G. and Tsang, W. W.
2330   * "A Simple Method for Generating Gamma Variables"
2331   * ACM Transactions on Mathematical Software, 26, 3, 363-372, 2000.
2332   */
2333  template<typename _RealType>
2334    template<typename _UniformRandomNumberGenerator>
2335      typename gamma_distribution<_RealType>::result_type
2336      gamma_distribution<_RealType>::
2337      operator()(_UniformRandomNumberGenerator& __urng,
2338  const param_type& __param)
2339      {
2340 __detail::_Adaptor<_UniformRandomNumberGenerator, result_type>
2341   __aurng(__urng);
2342
2343 result_type __u, __v, __n;
2344 const result_type __a1 = (__param._M_malpha
2345   - _RealType(1.0) / _RealType(3.0));
2346
2347 do
2348   {
2349     do
2350       {
2351 __n = _M_nd(__urng);
2352 __v = result_type(1.0) + __param._M_a2 * __n; 
2353       }
2354     while (__v <= 0.0);
2355
2356     __v = __v * __v * __v;
2357     __u = __aurng();
2358   }
2359 while (__u > result_type(1.0) - 0.0331 * __n * __n * __n * __n
2360        && (std::log(__u) > (0.5 * __n * __n + __a1
2361     * (1.0 - __v + std::log(__v)))));
2362
2363 if (__param.alpha() == __param._M_malpha)
2364   return __a1 * __v * __param.beta();
2365 else
2366   {
2367     do
2368       __u = __aurng();
2369     while (__u == 0.0);
2370     
2371     return (std::pow(__u, result_type(1.0) / __param.alpha())
2372     * __a1 * __v * __param.beta());
2373   }
2374      }
2375
2376  template<typename _RealType>
2377    template<typename _ForwardIterator,
2378      typename _UniformRandomNumberGenerator>
2379      void
2380      gamma_distribution<_RealType>::
2381      __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
2382       _UniformRandomNumberGenerator& __urng,
2383       const param_type& __param)
2384      {
2385 __glibcxx_function_requires(_ForwardIteratorConcept<_ForwardIterator>)
2386 __detail::_Adaptor<_UniformRandomNumberGenerator, result_type>
2387   __aurng(__urng);
2388
2389 result_type __u, __v, __n;
2390 const result_type __a1 = (__param._M_malpha
2391   - _RealType(1.0) / _RealType(3.0));
2392
2393 if (__param.alpha() == __param._M_malpha)
2394   while (__f != __t)
2395     {
2396       do
2397 {
2398   do
2399     {
2400       __n = _M_nd(__urng);
2401       __v = result_type(1.0) + __param._M_a2 * __n;
2402     }
2403   while (__v <= 0.0);
2404
2405   __v = __v * __v * __v;
2406   __u = __aurng();
2407 }
2408       while (__u > result_type(1.0) - 0.0331 * __n * __n * __n * __n
2409      && (std::log(__u) > (0.5 * __n * __n + __a1
2410   * (1.0 - __v + std::log(__v)))));
2411
2412       *__f++ = __a1 * __v * __param.beta();
2413     }
2414 else
2415   while (__f != __t)
2416     {
2417       do
2418 {
2419   do
2420     {
2421       __n = _M_nd(__urng);
2422       __v = result_type(1.0) + __param._M_a2 * __n;
2423     }
2424   while (__v <= 0.0);
2425
2426   __v = __v * __v * __v;
2427   __u = __aurng();
2428 }
2429       while (__u > result_type(1.0) - 0.0331 * __n * __n * __n * __n
2430      && (std::log(__u) > (0.5 * __n * __n + __a1
2431   * (1.0 - __v + std::log(__v)))));
2432
2433       do
2434 __u = __aurng();
2435       while (__u == 0.0);
2436
2437       *__f++ = (std::pow(__u, result_type(1.0) / __param.alpha())
2438 * __a1 * __v * __param.beta());
2439     }
2440      }
2441
2442  template<typename _RealType, typename _CharT, typename _Traits>
2443    std::basic_ostream<_CharT, _Traits>&
2444    operator<<(std::basic_ostream<_CharT, _Traits>& __os,
2445        const gamma_distribution<_RealType>& __x)
2446    {
2447      typedef std::basic_ostream<_CharT, _Traits>  __ostream_type;
2448      typedef typename __ostream_type::ios_base    __ios_base;
2449
2450      const typename __ios_base::fmtflags __flags = __os.flags();
2451      const _CharT __fill = __os.fill();
2452      const std::streamsize __precision = __os.precision();
2453      const _CharT __space = __os.widen(' ');
2454      __os.flags(__ios_base::scientific | __ios_base::left);
2455      __os.fill(__space);
2456      __os.precision(std::numeric_limits<_RealType>::max_digits10);
2457
2458      __os << __x.alpha() << __space << __x.beta()
2459    << __space << __x._M_nd;
2460
2461      __os.flags(__flags);
2462      __os.fill(__fill);
2463      __os.precision(__precision);
2464      return __os;
2465    }
2466
2467  template<typename _RealType, typename _CharT, typename _Traits>
2468    std::basic_istream<_CharT, _Traits>&
2469    operator>>(std::basic_istream<_CharT, _Traits>& __is,
2470        gamma_distribution<_RealType>& __x)
2471    {
2472      typedef std::basic_istream<_CharT, _Traits>  __istream_type;
2473      typedef typename __istream_type::ios_base    __ios_base;
2474
2475      const typename __ios_base::fmtflags __flags = __is.flags();
2476      __is.flags(__ios_base::dec | __ios_base::skipws);
2477
2478      _RealType __alpha_val, __beta_val;
2479      __is >> __alpha_val >> __beta_val >> __x._M_nd;
2480      __x.param(typename gamma_distribution<_RealType>::
2481 param_type(__alpha_val, __beta_val));
2482
2483      __is.flags(__flags);
2484      return __is;
2485    }
2486
2487
2488  template<typename _RealType>
2489    template<typename _UniformRandomNumberGenerator>
2490      typename weibull_distribution<_RealType>::result_type
2491      weibull_distribution<_RealType>::
2492      operator()(_UniformRandomNumberGenerator& __urng,
2493  const param_type& __p)
2494      {
2495 __detail::_Adaptor<_UniformRandomNumberGenerator, result_type>
2496   __aurng(__urng);
2497 return __p.b() * std::pow(-std::log(result_type(1) - __aurng()),
2498   result_type(1) / __p.a());
2499      }
2500
2501  template<typename _RealType>
2502    template<typename _ForwardIterator,
2503      typename _UniformRandomNumberGenerator>
2504      void
2505      weibull_distribution<_RealType>::
2506      __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
2507       _UniformRandomNumberGenerator& __urng,
2508       const param_type& __p)
2509      {
2510 __glibcxx_function_requires(_ForwardIteratorConcept<_ForwardIterator>)
2511 __detail::_Adaptor<_UniformRandomNumberGenerator, result_type>
2512   __aurng(__urng);
2513 auto __inv_a = result_type(1) / __p.a();
2514
2515 while (__f != __t)
2516   *__f++ = __p.b() * std::pow(-std::log(result_type(1) - __aurng()),
2517       __inv_a);
2518      }
2519
2520  template<typename _RealType, typename _CharT, typename _Traits>
2521    std::basic_ostream<_CharT, _Traits>&
2522    operator<<(std::basic_ostream<_CharT, _Traits>& __os,
2523        const weibull_distribution<_RealType>& __x)
2524    {
2525      typedef std::basic_ostream<_CharT, _Traits>  __ostream_type;
2526      typedef typename __ostream_type::ios_base    __ios_base;
2527
2528      const typename __ios_base::fmtflags __flags = __os.flags();
2529      const _CharT __fill = __os.fill();
2530      const std::streamsize __precision = __os.precision();
2531      const _CharT __space = __os.widen(' ');
2532      __os.flags(__ios_base::scientific | __ios_base::left);
2533      __os.fill(__space);
2534      __os.precision(std::numeric_limits<_RealType>::max_digits10);
2535
2536      __os << __x.a() << __space << __x.b();
2537
2538      __os.flags(__flags);
2539      __os.fill(__fill);
2540      __os.precision(__precision);
2541      return __os;
2542    }
2543
2544  template<typename _RealType, typename _CharT, typename _Traits>
2545    std::basic_istream<_CharT, _Traits>&
2546    operator>>(std::basic_istream<_CharT, _Traits>& __is,
2547        weibull_distribution<_RealType>& __x)
2548    {
2549      typedef std::basic_istream<_CharT, _Traits>  __istream_type;
2550      typedef typename __istream_type::ios_base    __ios_base;
2551
2552      const typename __ios_base::fmtflags __flags = __is.flags();
2553      __is.flags(__ios_base::dec | __ios_base::skipws);
2554
2555      _RealType __a, __b;
2556      __is >> __a >> __b;
2557      __x.param(typename weibull_distribution<_RealType>::
2558 param_type(__a, __b));
2559
2560      __is.flags(__flags);
2561      return __is;
2562    }
2563
2564
2565  template<typename _RealType>
2566    template<typename _UniformRandomNumberGenerator>
2567      typename extreme_value_distribution<_RealType>::result_type
2568      extreme_value_distribution<_RealType>::
2569      operator()(_UniformRandomNumberGenerator& __urng,
2570  const param_type& __p)
2571      {
2572 __detail::_Adaptor<_UniformRandomNumberGenerator, result_type>
2573   __aurng(__urng);
2574 return __p.a() - __p.b() * std::log(-std::log(result_type(1)
2575       - __aurng()));
2576      }
2577
2578  template<typename _RealType>
2579    template<typename _ForwardIterator,
2580      typename _UniformRandomNumberGenerator>
2581      void
2582      extreme_value_distribution<_RealType>::
2583      __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
2584       _UniformRandomNumberGenerator& __urng,
2585       const param_type& __p)
2586      {
2587 __glibcxx_function_requires(_ForwardIteratorConcept<_ForwardIterator>)
2588 __detail::_Adaptor<_UniformRandomNumberGenerator, result_type>
2589   __aurng(__urng);
2590
2591 while (__f != __t)
2592   *__f++ = __p.a() - __p.b() * std::log(-std::log(result_type(1)
2593   - __aurng()));
2594      }
2595
2596  template<typename _RealType, typename _CharT, typename _Traits>
2597    std::basic_ostream<_CharT, _Traits>&
2598    operator<<(std::basic_ostream<_CharT, _Traits>& __os,
2599        const extreme_value_distribution<_RealType>& __x)
2600    {
2601      typedef std::basic_ostream<_CharT, _Traits>  __ostream_type;
2602      typedef typename __ostream_type::ios_base    __ios_base;
2603
2604      const typename __ios_base::fmtflags __flags = __os.flags();
2605      const _CharT __fill = __os.fill();
2606      const std::streamsize __precision = __os.precision();
2607      const _CharT __space = __os.widen(' ');
2608      __os.flags(__ios_base::scientific | __ios_base::left);
2609      __os.fill(__space);
2610      __os.precision(std::numeric_limits<_RealType>::max_digits10);
2611
2612      __os << __x.a() << __space << __x.b();
2613
2614      __os.flags(__flags);
2615      __os.fill(__fill);
2616      __os.precision(__precision);
2617      return __os;
2618    }
2619
2620  template<typename _RealType, typename _CharT, typename _Traits>
2621    std::basic_istream<_CharT, _Traits>&
2622    operator>>(std::basic_istream<_CharT, _Traits>& __is,
2623        extreme_value_distribution<_RealType>& __x)
2624    {
2625      typedef std::basic_istream<_CharT, _Traits>  __istream_type;
2626      typedef typename __istream_type::ios_base    __ios_base;
2627
2628      const typename __ios_base::fmtflags __flags = __is.flags();
2629      __is.flags(__ios_base::dec | __ios_base::skipws);
2630
2631      _RealType __a, __b;
2632      __is >> __a >> __b;
2633      __x.param(typename extreme_value_distribution<_RealType>::
2634 param_type(__a, __b));
2635
2636      __is.flags(__flags);
2637      return __is;
2638    }
2639
2640
2641  template<typename _IntType>
2642    void
2643    discrete_distribution<_IntType>::param_type::
2644    _M_initialize()
2645    {
2646      if (_M_prob.size() < 2)
2647 {
2648   _M_prob.clear();
2649   return;
2650 }
2651
2652      const double __sum = std::accumulate(_M_prob.begin(),
2653    _M_prob.end(), 0.0);
2654      // Now normalize the probabilites.
2655      __detail::__normalize(_M_prob.begin(), _M_prob.end(), _M_prob.begin(),
2656     __sum);
2657      // Accumulate partial sums.
2658      _M_cp.reserve(_M_prob.size());
2659      std::partial_sum(_M_prob.begin(), _M_prob.end(),
2660        std::back_inserter(_M_cp));
2661      // Make sure the last cumulative probability is one.
2662      _M_cp[_M_cp.size() - 1] = 1.0;
2663    }
2664
2665  template<typename _IntType>
2666    template<typename _Func>
2667      discrete_distribution<_IntType>::param_type::
2668      param_type(size_t __nw, double __xmin, double __xmax, _Func __fw)
2669      : _M_prob(), _M_cp()
2670      {
2671 const size_t __n = __nw == 0 ? 1 : __nw;
2672 const double __delta = (__xmax - __xmin) / __n;
2673
2674 _M_prob.reserve(__n);
2675 for (size_t __k = 0; __k < __nw; ++__k)
2676   _M_prob.push_back(__fw(__xmin + __k * __delta + 0.5 * __delta));
2677
2678 _M_initialize();
2679      }
2680
2681  template<typename _IntType>
2682    template<typename _UniformRandomNumberGenerator>
2683      typename discrete_distribution<_IntType>::result_type
2684      discrete_distribution<_IntType>::
2685      operator()(_UniformRandomNumberGenerator& __urng,
2686  const param_type& __param)
2687      {
2688 if (__param._M_cp.empty())
2689   return result_type(0);
2690
2691 __detail::_Adaptor<_UniformRandomNumberGenerator, double>
2692   __aurng(__urng);
2693
2694 const double __p = __aurng();
2695 auto __pos = std::lower_bound(__param._M_cp.begin(),
2696       __param._M_cp.end(), __p);
2697
2698 return __pos - __param._M_cp.begin();
2699      }
2700
2701  template<typename _IntType>
2702    template<typename _ForwardIterator,
2703      typename _UniformRandomNumberGenerator>
2704      void
2705      discrete_distribution<_IntType>::
2706      __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
2707       _UniformRandomNumberGenerator& __urng,
2708       const param_type& __param)
2709      {
2710 __glibcxx_function_requires(_ForwardIteratorConcept<_ForwardIterator>)
2711
2712 if (__param._M_cp.empty())
2713   {
2714     while (__f != __t)
2715       *__f++ = result_type(0);
2716     return;
2717   }
2718
2719 __detail::_Adaptor<_UniformRandomNumberGenerator, double>
2720   __aurng(__urng);
2721
2722 while (__f != __t)
2723   {
2724     const double __p = __aurng();
2725     auto __pos = std::lower_bound(__param._M_cp.begin(),
2726   __param._M_cp.end(), __p);
2727
2728     *__f++ = __pos - __param._M_cp.begin();
2729   }
2730      }
2731
2732  template<typename _IntType, typename _CharT, typename _Traits>
2733    std::basic_ostream<_CharT, _Traits>&
2734    operator<<(std::basic_ostream<_CharT, _Traits>& __os,
2735        const discrete_distribution<_IntType>& __x)
2736    {
2737      typedef std::basic_ostream<_CharT, _Traits>  __ostream_type;
2738      typedef typename __ostream_type::ios_base    __ios_base;
2739
2740      const typename __ios_base::fmtflags __flags = __os.flags();
2741      const _CharT __fill = __os.fill();
2742      const std::streamsize __precision = __os.precision();
2743      const _CharT __space = __os.widen(' ');
2744      __os.flags(__ios_base::scientific | __ios_base::left);
2745      __os.fill(__space);
2746      __os.precision(std::numeric_limits<double>::max_digits10);
2747
2748      std::vector<double> __prob = __x.probabilities();
2749      __os << __prob.size();
2750      for (auto __dit = __prob.begin(); __dit != __prob.end(); ++__dit)
2751 __os << __space << *__dit;
2752
2753      __os.flags(__flags);
2754      __os.fill(__fill);
2755      __os.precision(__precision);
2756      return __os;
2757    }
2758
2759  template<typename _IntType, typename _CharT, typename _Traits>
2760    std::basic_istream<_CharT, _Traits>&
2761    operator>>(std::basic_istream<_CharT, _Traits>& __is,
2762        discrete_distribution<_IntType>& __x)
2763    {
2764      typedef std::basic_istream<_CharT, _Traits>  __istream_type;
2765      typedef typename __istream_type::ios_base    __ios_base;
2766
2767      const typename __ios_base::fmtflags __flags = __is.flags();
2768      __is.flags(__ios_base::dec | __ios_base::skipws);
2769
2770      size_t __n;
2771      __is >> __n;
2772
2773      std::vector<double> __prob_vec;
2774      __prob_vec.reserve(__n);
2775      for (; __n != 0; --__n)
2776 {
2777   double __prob;
2778   __is >> __prob;
2779   __prob_vec.push_back(__prob);
2780 }
2781
2782      __x.param(typename discrete_distribution<_IntType>::
2783 param_type(__prob_vec.begin(), __prob_vec.end()));
2784
2785      __is.flags(__flags);
2786      return __is;
2787    }
2788
2789
2790  template<typename _RealType>
2791    void
2792    piecewise_constant_distribution<_RealType>::param_type::
2793    _M_initialize()
2794    {
2795      if (_M_int.size() < 2
2796   || (_M_int.size() == 2
2797       && _M_int[0] == _RealType(0)
2798       && _M_int[1] == _RealType(1)))
2799 {
2800   _M_int.clear();
2801   _M_den.clear();
2802   return;
2803 }
2804
2805      const double __sum = std::accumulate(_M_den.begin(),
2806    _M_den.end(), 0.0);
2807
2808      __detail::__normalize(_M_den.begin(), _M_den.end(), _M_den.begin(),
2809     __sum);
2810
2811      _M_cp.reserve(_M_den.size());
2812      std::partial_sum(_M_den.begin(), _M_den.end(),
2813        std::back_inserter(_M_cp));
2814
2815      // Make sure the last cumulative probability is one.
2816      _M_cp[_M_cp.size() - 1] = 1.0;
2817
2818      for (size_t __k = 0; __k < _M_den.size(); ++__k)
2819 _M_den[__k] /= _M_int[__k + 1] - _M_int[__k];
2820    }
2821
2822  template<typename _RealType>
2823    template<typename _InputIteratorB, typename _InputIteratorW>
2824      piecewise_constant_distribution<_RealType>::param_type::
2825      param_type(_InputIteratorB __bbegin,
2826  _InputIteratorB __bend,
2827  _InputIteratorW __wbegin)
2828      : _M_int(), _M_den(), _M_cp()
2829      {
2830 if (__bbegin != __bend)
2831   {
2832     for (;;)
2833       {
2834 _M_int.push_back(*__bbegin);
2835 ++__bbegin;
2836 if (__bbegin == __bend)
2837   break;
2838
2839 _M_den.push_back(*__wbegin);
2840 ++__wbegin;
2841       }
2842   }
2843
2844 _M_initialize();
2845      }
2846
2847  template<typename _RealType>
2848    template<typename _Func>
2849      piecewise_constant_distribution<_RealType>::param_type::
2850      param_type(initializer_list<_RealType> __bl, _Func __fw)
2851      : _M_int(), _M_den(), _M_cp()
2852      {
2853 _M_int.reserve(__bl.size());
2854 for (auto __biter = __bl.begin(); __biter != __bl.end(); ++__biter)
2855   _M_int.push_back(*__biter);
2856
2857 _M_den.reserve(_M_int.size() - 1);
2858 for (size_t __k = 0; __k < _M_int.size() - 1; ++__k)
2859   _M_den.push_back(__fw(0.5 * (_M_int[__k + 1] + _M_int[__k])));
2860
2861 _M_initialize();
2862      }
2863
2864  template<typename _RealType>
2865    template<typename _Func>
2866      piecewise_constant_distribution<_RealType>::param_type::
2867      param_type(size_t __nw, _RealType __xmin, _RealType __xmax, _Func __fw)
2868      : _M_int(), _M_den(), _M_cp()
2869      {
2870 const size_t __n = __nw == 0 ? 1 : __nw;
2871 const _RealType __delta = (__xmax - __xmin) / __n;
2872
2873 _M_int.reserve(__n + 1);
2874 for (size_t __k = 0; __k <= __nw; ++__k)
2875   _M_int.push_back(__xmin + __k * __delta);
2876
2877 _M_den.reserve(__n);
2878 for (size_t __k = 0; __k < __nw; ++__k)
2879   _M_den.push_back(__fw(_M_int[__k] + 0.5 * __delta));
2880
2881 _M_initialize();
2882      }
2883
2884  template<typename _RealType>
2885    template<typename _UniformRandomNumberGenerator>
2886      typename piecewise_constant_distribution<_RealType>::result_type
2887      piecewise_constant_distribution<_RealType>::
2888      operator()(_UniformRandomNumberGenerator& __urng,
2889  const param_type& __param)
2890      {
2891 __detail::_Adaptor<_UniformRandomNumberGenerator, double>
2892   __aurng(__urng);
2893
2894 const double __p = __aurng();
2895 if (__param._M_cp.empty())
2896   return __p;
2897
2898 auto __pos = std::lower_bound(__param._M_cp.begin(),
2899       __param._M_cp.end(), __p);
2900 const size_t __i = __pos - __param._M_cp.begin();
2901
2902 const double __pref = __i > 0 ? __param._M_cp[__i - 1] : 0.0;
2903
2904 return __param._M_int[__i] + (__p - __pref) / __param._M_den[__i];
2905      }
2906
2907  template<typename _RealType>
2908    template<typename _ForwardIterator,
2909      typename _UniformRandomNumberGenerator>
2910      void
2911      piecewise_constant_distribution<_RealType>::
2912      __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
2913       _UniformRandomNumberGenerator& __urng,
2914       const param_type& __param)
2915      {
2916 __glibcxx_function_requires(_ForwardIteratorConcept<_ForwardIterator>)
2917 __detail::_Adaptor<_UniformRandomNumberGenerator, double>
2918   __aurng(__urng);
2919
2920 if (__param._M_cp.empty())
2921   {
2922     while (__f != __t)
2923       *__f++ = __aurng();
2924     return;
2925   }
2926
2927 while (__f != __t)
2928   {
2929     const double __p = __aurng();
2930
2931     auto __pos = std::lower_bound(__param._M_cp.begin(),
2932   __param._M_cp.end(), __p);
2933     const size_t __i = __pos - __param._M_cp.begin();
2934
2935     const double __pref = __i > 0 ? __param._M_cp[__i - 1] : 0.0;
2936
2937     *__f++ = (__param._M_int[__i]
2938       + (__p - __pref) / __param._M_den[__i]);
2939   }
2940      }
2941
2942  template<typename _RealType, typename _CharT, typename _Traits>
2943    std::basic_ostream<_CharT, _Traits>&
2944    operator<<(std::basic_ostream<_CharT, _Traits>& __os,
2945        const piecewise_constant_distribution<_RealType>& __x)
2946    {
2947      typedef std::basic_ostream<_CharT, _Traits>  __ostream_type;
2948      typedef typename __ostream_type::ios_base    __ios_base;
2949
2950      const typename __ios_base::fmtflags __flags = __os.flags();
2951      const _CharT __fill = __os.fill();
2952      const std::streamsize __precision = __os.precision();
2953      const _CharT __space = __os.widen(' ');
2954      __os.flags(__ios_base::scientific | __ios_base::left);
2955      __os.fill(__space);
2956      __os.precision(std::numeric_limits<_RealType>::max_digits10);
2957
2958      std::vector<_RealType> __int = __x.intervals();
2959      __os << __int.size() - 1;
2960
2961      for (auto __xit = __int.begin(); __xit != __int.end(); ++__xit)
2962 __os << __space << *__xit;
2963
2964      std::vector<double> __den = __x.densities();
2965      for (auto __dit = __den.begin(); __dit != __den.end(); ++__dit)
2966 __os << __space << *__dit;
2967
2968      __os.flags(__flags);
2969      __os.fill(__fill);
2970      __os.precision(__precision);
2971      return __os;
2972    }
2973
2974  template<typename _RealType, typename _CharT, typename _Traits>
2975    std::basic_istream<_CharT, _Traits>&
2976    operator>>(std::basic_istream<_CharT, _Traits>& __is,
2977        piecewise_constant_distribution<_RealType>& __x)
2978    {
2979      typedef std::basic_istream<_CharT, _Traits>  __istream_type;
2980      typedef typename __istream_type::ios_base    __ios_base;
2981
2982      const typename __ios_base::fmtflags __flags = __is.flags();
2983      __is.flags(__ios_base::dec | __ios_base::skipws);
2984
2985      size_t __n;
2986      __is >> __n;
2987
2988      std::vector<_RealType> __int_vec;
2989      __int_vec.reserve(__n + 1);
2990      for (size_t __i = 0; __i <= __n; ++__i)
2991 {
2992   _RealType __int;
2993   __is >> __int;
2994   __int_vec.push_back(__int);
2995 }
2996
2997      std::vector<double> __den_vec;
2998      __den_vec.reserve(__n);
2999      for (size_t __i = 0; __i < __n; ++__i)
3000 {
3001   double __den;
3002   __is >> __den;
3003   __den_vec.push_back(__den);
3004 }
3005
3006      __x.param(typename piecewise_constant_distribution<_RealType>::
3007   param_type(__int_vec.begin(), __int_vec.end(), __den_vec.begin()));
3008
3009      __is.flags(__flags);
3010      return __is;
3011    }
3012
3013
3014  template<typename _RealType>
3015    void
3016    piecewise_linear_distribution<_RealType>::param_type::
3017    _M_initialize()
3018    {
3019      if (_M_int.size() < 2
3020   || (_M_int.size() == 2
3021       && _M_int[0] == _RealType(0)
3022       && _M_int[1] == _RealType(1)
3023       && _M_den[0] == _M_den[1]))
3024 {
3025   _M_int.clear();
3026   _M_den.clear();
3027   return;
3028 }
3029
3030      double __sum = 0.0;
3031      _M_cp.reserve(_M_int.size() - 1);
3032      _M_m.reserve(_M_int.size() - 1);
3033      for (size_t __k = 0; __k < _M_int.size() - 1; ++__k)
3034 {
3035   const _RealType __delta = _M_int[__k + 1] - _M_int[__k];
3036   __sum += 0.5 * (_M_den[__k + 1] + _M_den[__k]) * __delta;
3037   _M_cp.push_back(__sum);
3038   _M_m.push_back((_M_den[__k + 1] - _M_den[__k]) / __delta);
3039 }
3040
3041      //  Now normalize the densities...
3042      __detail::__normalize(_M_den.begin(), _M_den.end(), _M_den.begin(),
3043     __sum);
3044      //  ... and partial sums... 
3045      __detail::__normalize(_M_cp.begin(), _M_cp.end(), _M_cp.begin(), __sum);
3046      //  ... and slopes.
3047      __detail::__normalize(_M_m.begin(), _M_m.end(), _M_m.begin(), __sum);
3048
3049      //  Make sure the last cumulative probablility is one.
3050      _M_cp[_M_cp.size() - 1] = 1.0;
3051     }
3052
3053  template<typename _RealType>
3054    template<typename _InputIteratorB, typename _InputIteratorW>
3055      piecewise_linear_distribution<_RealType>::param_type::
3056      param_type(_InputIteratorB __bbegin,
3057  _InputIteratorB __bend,
3058  _InputIteratorW __wbegin)
3059      : _M_int(), _M_den(), _M_cp(), _M_m()
3060      {
3061 for (; __bbegin != __bend; ++__bbegin, ++__wbegin)
3062   {
3063     _M_int.push_back(*__bbegin);
3064     _M_den.push_back(*__wbegin);
3065   }
3066
3067 _M_initialize();
3068      }
3069
3070  template<typename _RealType>
3071    template<typename _Func>
3072      piecewise_linear_distribution<_RealType>::param_type::
3073      param_type(initializer_list<_RealType> __bl, _Func __fw)
3074      : _M_int(), _M_den(), _M_cp(), _M_m()
3075      {
3076 _M_int.reserve(__bl.size());
3077 _M_den.reserve(__bl.size());
3078 for (auto __biter = __bl.begin(); __biter != __bl.end(); ++__biter)
3079   {
3080     _M_int.push_back(*__biter);
3081     _M_den.push_back(__fw(*__biter));
3082   }
3083
3084 _M_initialize();
3085      }
3086
3087  template<typename _RealType>
3088    template<typename _Func>
3089      piecewise_linear_distribution<_RealType>::param_type::
3090      param_type(size_t __nw, _RealType __xmin, _RealType __xmax, _Func __fw)
3091      : _M_int(), _M_den(), _M_cp(), _M_m()
3092      {
3093 const size_t __n = __nw == 0 ? 1 : __nw;
3094 const _RealType __delta = (__xmax - __xmin) / __n;
3095
3096 _M_int.reserve(__n + 1);
3097 _M_den.reserve(__n + 1);
3098 for (size_t __k = 0; __k <= __nw; ++__k)
3099   {
3100     _M_int.push_back(__xmin + __k * __delta);
3101     _M_den.push_back(__fw(_M_int[__k] + __delta));
3102   }
3103
3104 _M_initialize();
3105      }
3106
3107  template<typename _RealType>
3108    template<typename _UniformRandomNumberGenerator>
3109      typename piecewise_linear_distribution<_RealType>::result_type
3110      piecewise_linear_distribution<_RealType>::
3111      operator()(_UniformRandomNumberGenerator& __urng,
3112  const param_type& __param)
3113      {
3114 __detail::_Adaptor<_UniformRandomNumberGenerator, double>
3115   __aurng(__urng);
3116
3117 const double __p = __aurng();
3118 if (__param._M_cp.empty())
3119   return __p;
3120
3121 auto __pos = std::lower_bound(__param._M_cp.begin(),
3122       __param._M_cp.end(), __p);
3123 const size_t __i = __pos - __param._M_cp.begin();
3124
3125 const double __pref = __i > 0 ? __param._M_cp[__i - 1] : 0.0;
3126
3127 const double __a = 0.5 * __param._M_m[__i];
3128 const double __b = __param._M_den[__i];
3129 const double __cm = __p - __pref;
3130
3131 _RealType __x = __param._M_int[__i];
3132 if (__a == 0)
3133   __x += __cm / __b;
3134 else
3135   {
3136     const double __d = __b * __b + 4.0 * __a * __cm;
3137     __x += 0.5 * (std::sqrt(__d) - __b) / __a;
3138          }
3139
3140        return __x;
3141      }
3142
3143  template<typename _RealType>
3144    template<typename _ForwardIterator,
3145      typename _UniformRandomNumberGenerator>
3146      void
3147      piecewise_linear_distribution<_RealType>::
3148      __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
3149       _UniformRandomNumberGenerator& __urng,
3150       const param_type& __param)
3151      {
3152 __glibcxx_function_requires(_ForwardIteratorConcept<_ForwardIterator>)
3153 // We could duplicate everything from operator()...
3154 while (__f != __t)
3155   *__f++ = this->operator()(__urng, __param);
3156      }
3157
3158  template<typename _RealType, typename _CharT, typename _Traits>
3159    std::basic_ostream<_CharT, _Traits>&
3160    operator<<(std::basic_ostream<_CharT, _Traits>& __os,
3161        const piecewise_linear_distribution<_RealType>& __x)
3162    {
3163      typedef std::basic_ostream<_CharT, _Traits>  __ostream_type;
3164      typedef typename __ostream_type::ios_base    __ios_base;
3165
3166      const typename __ios_base::fmtflags __flags = __os.flags();
3167      const _CharT __fill = __os.fill();
3168      const std::streamsize __precision = __os.precision();
3169      const _CharT __space = __os.widen(' ');
3170      __os.flags(__ios_base::scientific | __ios_base::left);
3171      __os.fill(__space);
3172      __os.precision(std::numeric_limits<_RealType>::max_digits10);
3173
3174      std::vector<_RealType> __int = __x.intervals();
3175      __os << __int.size() - 1;
3176
3177      for (auto __xit = __int.begin(); __xit != __int.end(); ++__xit)
3178 __os << __space << *__xit;
3179
3180      std::vector<double> __den = __x.densities();
3181      for (auto __dit = __den.begin(); __dit != __den.end(); ++__dit)
3182 __os << __space << *__dit;
3183
3184      __os.flags(__flags);
3185      __os.fill(__fill);
3186      __os.precision(__precision);
3187      return __os;
3188    }
3189
3190  template<typename _RealType, typename _CharT, typename _Traits>
3191    std::basic_istream<_CharT, _Traits>&
3192    operator>>(std::basic_istream<_CharT, _Traits>& __is,
3193        piecewise_linear_distribution<_RealType>& __x)
3194    {
3195      typedef std::basic_istream<_CharT, _Traits>  __istream_type;
3196      typedef typename __istream_type::ios_base    __ios_base;
3197
3198      const typename __ios_base::fmtflags __flags = __is.flags();
3199      __is.flags(__ios_base::dec | __ios_base::skipws);
3200
3201      size_t __n;
3202      __is >> __n;
3203
3204      std::vector<_RealType> __int_vec;
3205      __int_vec.reserve(__n + 1);
3206      for (size_t __i = 0; __i <= __n; ++__i)
3207 {
3208   _RealType __int;
3209   __is >> __int;
3210   __int_vec.push_back(__int);
3211 }
3212
3213      std::vector<double> __den_vec;
3214      __den_vec.reserve(__n + 1);
3215      for (size_t __i = 0; __i <= __n; ++__i)
3216 {
3217   double __den;
3218   __is >> __den;
3219   __den_vec.push_back(__den);
3220 }
3221
3222      __x.param(typename piecewise_linear_distribution<_RealType>::
3223   param_type(__int_vec.begin(), __int_vec.end(), __den_vec.begin()));
3224
3225      __is.flags(__flags);
3226      return __is;
3227    }
3228
3229
3230  template<typename _IntType>
3231    seed_seq::seed_seq(std::initializer_list<_IntType> __il)
3232    {
3233      for (auto __iter = __il.begin(); __iter != __il.end(); ++__iter)
3234 _M_v.push_back(__detail::__mod<result_type,
3235        __detail::_Shift<result_type, 32>::__value>(*__iter));
3236    }
3237
3238  template<typename _InputIterator>
3239    seed_seq::seed_seq(_InputIterator __begin, _InputIterator __end)
3240    {
3241      for (_InputIterator __iter = __begin; __iter != __end; ++__iter)
3242 _M_v.push_back(__detail::__mod<result_type,
3243        __detail::_Shift<result_type, 32>::__value>(*__iter));
3244    }
3245
3246  template<typename _RandomAccessIterator>
3247    void
3248    seed_seq::generate(_RandomAccessIterator __begin,
3249        _RandomAccessIterator __end)
3250    {
3251      typedef typename iterator_traits<_RandomAccessIterator>::value_type
3252        _Type;
3253
3254      if (__begin == __end)
3255 return;
3256
3257      std::fill(__begin, __end, _Type(0x8b8b8b8bu));
3258
3259      const size_t __n = __end - __begin;
3260      const size_t __s = _M_v.size();
3261      const size_t __t = (__n >= 623) ? 11
3262        : (__n >=  68) ? 7
3263        : (__n >=  39) ? 5
3264        : (__n >=   7) ? 3
3265        : (__n - 1) / 2;
3266      const size_t __p = (__n - __t) / 2;
3267      const size_t __q = __p + __t;
3268      const size_t __m = std::max(size_t(__s + 1), __n);
3269
3270      for (size_t __k = 0; __k < __m; ++__k)
3271 {
3272   _Type __arg = (__begin[__k % __n]
3273  ^ __begin[(__k + __p) % __n]
3274  ^ __begin[(__k - 1) % __n]);
3275   _Type __r1 = __arg ^ (__arg >> 27);
3276   __r1 = __detail::__mod<_Type,
3277     __detail::_Shift<_Type, 32>::__value>(1664525u * __r1);
3278   _Type __r2 = __r1;
3279   if (__k == 0)
3280     __r2 += __s;
3281   else if (__k <= __s)
3282     __r2 += __k % __n + _M_v[__k - 1];
3283   else
3284     __r2 += __k % __n;
3285   __r2 = __detail::__mod<_Type,
3286            __detail::_Shift<_Type, 32>::__value>(__r2);
3287   __begin[(__k + __p) % __n] += __r1;
3288   __begin[(__k + __q) % __n] += __r2;
3289   __begin[__k % __n] = __r2;
3290 }
3291
3292      for (size_t __k = __m; __k < __m + __n; ++__k)
3293 {
3294   _Type __arg = (__begin[__k % __n]
3295  + __begin[(__k + __p) % __n]
3296  + __begin[(__k - 1) % __n]);
3297   _Type __r3 = __arg ^ (__arg >> 27);
3298   __r3 = __detail::__mod<_Type,
3299    __detail::_Shift<_Type, 32>::__value>(1566083941u * __r3);
3300   _Type __r4 = __r3 - __k % __n;
3301   __r4 = __detail::__mod<_Type,
3302            __detail::_Shift<_Type, 32>::__value>(__r4);
3303   __begin[(__k + __p) % __n] ^= __r3;
3304   __begin[(__k + __q) % __n] ^= __r4;
3305   __begin[__k % __n] = __r4;
3306 }
3307    }
3308
3309  template<typename _RealType, size_t __bits,
3310    typename _UniformRandomNumberGenerator>
3311    _RealType
3312    generate_canonical(_UniformRandomNumberGenerator& __urng)
3313    {
3314      static_assert(std::is_floating_point<_RealType>::value,
3315     "template argument must be a floating point type");
3316
3317      const size_t __b
3318 = std::min(static_cast<size_t>(std::numeric_limits<_RealType>::digits),
3319                   __bits);
3320      const long double __r = static_cast<long double>(__urng.max())
3321     - static_cast<long double>(__urng.min()) + 1.0L;
3322      const size_t __log2r = std::log(__r) / std::log(2.0L);
3323      const size_t __m = std::max<size_t>(1UL,
3324   (__b + __log2r - 1UL) / __log2r);
3325      _RealType __ret;
3326      _RealType __sum = _RealType(0);
3327      _RealType __tmp = _RealType(1);
3328      for (size_t __k = __m; __k != 0; --__k)
3329 {
3330   __sum += _RealType(__urng() - __urng.min()) * __tmp;
3331   __tmp *= __r;
3332 }
3333      __ret = __sum / __tmp;
3334      if (__builtin_expect(__ret >= _RealType(1), 0))
3335 {
3336#if _GLIBCXX_USE_C99_MATH_TR1
3337   __ret = std::nextafter(_RealType(1), _RealType(0));
3338#else
3339   __ret = _RealType(1)
3340     - std::numeric_limits<_RealType>::epsilon() / _RealType(2);
3341#endif
3342 }
3343      return __ret;
3344    }
3345
3346_GLIBCXX_END_NAMESPACE_VERSION
3347} // namespace
3348
3349#endif
3350
std::__detail::_Mod::__calc
std::linear_congruential_engine::multiplier
std::linear_congruential_engine::increment
std::linear_congruential_engine::modulus
std::linear_congruential_engine::default_seed
std::linear_congruential_engine::seed
std::linear_congruential_engine::seed
std::mersenne_twister_engine::word_size
std::mersenne_twister_engine::state_size
std::mersenne_twister_engine::shift_size
std::mersenne_twister_engine::mask_bits
std::mersenne_twister_engine::xor_mask
std::mersenne_twister_engine::tempering_u
std::mersenne_twister_engine::tempering_d
std::mersenne_twister_engine::tempering_s
std::mersenne_twister_engine::tempering_b
std::mersenne_twister_engine::tempering_t
std::mersenne_twister_engine::tempering_c
std::mersenne_twister_engine::tempering_l
std::mersenne_twister_engine::initialization_multiplier
std::mersenne_twister_engine::default_seed
std::mersenne_twister_engine::seed
std::mersenne_twister_engine::seed
std::mersenne_twister_engine::_M_gen_rand
std::mersenne_twister_engine::discard
std::subtract_with_carry_engine::word_size
std::subtract_with_carry_engine::short_lag
std::subtract_with_carry_engine::long_lag
std::subtract_with_carry_engine::default_seed
std::subtract_with_carry_engine::seed
std::subtract_with_carry_engine::seed
std::discard_block_engine::block_size
std::discard_block_engine::used_block
std::shuffle_order_engine::table_size
std::uniform_real_distribution::__generate_impl
std::bernoulli_distribution::__generate_impl
std::geometric_distribution::__generate_impl
std::negative_binomial_distribution::__generate_impl
std::negative_binomial_distribution::__generate_impl
std::poisson_distribution::param_type::_M_initialize
std::poisson_distribution::__generate_impl
std::binomial_distribution::param_type::_M_initialize
std::binomial_distribution::_M_waiting
std::binomial_distribution::__generate_impl
std::exponential_distribution::__generate_impl
std::normal_distribution::__generate_impl
std::lognormal_distribution::__generate_impl
std::chi_squared_distribution::__generate_impl
std::chi_squared_distribution::__generate_impl
std::cauchy_distribution::__generate_impl
std::fisher_f_distribution::__generate_impl
std::fisher_f_distribution::__generate_impl
std::student_t_distribution::__generate_impl
std::student_t_distribution::__generate_impl
std::gamma_distribution::param_type::_M_initialize
std::gamma_distribution::__generate_impl
std::weibull_distribution::__generate_impl
std::extreme_value_distribution::__generate_impl
std::discrete_distribution::param_type::_M_initialize
std::discrete_distribution::__generate_impl
std::piecewise_constant_distribution::param_type::_M_initialize
std::piecewise_constant_distribution::__generate_impl
std::piecewise_linear_distribution::param_type::_M_initialize
std::piecewise_linear_distribution::__generate_impl
std::seed_seq::generate
CodeThinker, Enjoy your reading code