libstdc++-v3/include/bits/random.h

   1 // random number generation -*- C++ -*-
   2
   3 // Copyright (C) 2009, 2010 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 /**
  26  * @file bits/random.h
  27  *  This is an internal header file, included by other library headers.
  28  *  You should not attempt to use it directly.
  29  */
  30
  31 #include <vector>
  32
  33 namespace std
  34 {
  35   // [26.4] Random number generation
  36
  37   /**
  38    * @defgroup random Random Number Generation
  39    * @ingroup numerics
  40    *
  41    * A facility for generating random numbers on selected distributions.
  42    * @{
  43    */
  44
  45   /**
  46    * @brief A function template for converting the output of a (integral)
  47    * uniform random number generator to a floatng point result in the range
  48    * [0-1).
  49    */
  50   template<typename _RealType, size_t __bits,
  51            typename _UniformRandomNumberGenerator>
  52     _RealType
  53     generate_canonical(_UniformRandomNumberGenerator& __g);
  54
  55   /*
  56    * Implementation-space details.
  57    */
  58   namespace __detail
  59   {
  60     template<typename _UIntType, size_t __w,
  61              bool = __w < static_cast<size_t>
  62                           (std::numeric_limits<_UIntType>::digits)>
  63       struct _Shift
  64       { static const _UIntType __value = 0; };
  65
  66     template<typename _UIntType, size_t __w>
  67       struct _Shift<_UIntType, __w, true>
  68       { static const _UIntType __value = _UIntType(1) << __w; };
  69
  70     template<typename _Tp, _Tp __m, _Tp __a, _Tp __c, bool>
  71       struct _Mod;
  72
  73     // Dispatch based on modulus value to prevent divide-by-zero compile-time
  74     // errors when m == 0.
  75     template<typename _Tp, _Tp __m, _Tp __a = 1, _Tp __c = 0>
  76       inline _Tp
  77       __mod(_Tp __x)
  78       { return _Mod<_Tp, __m, __a, __c, __m == 0>::__calc(__x); }
  79
  80     /*
  81      * An adaptor class for converting the output of any Generator into
  82      * the input for a specific Distribution.
  83      */
  84     template<typename _Engine, typename _DInputType>
  85       struct _Adaptor
  86       {
  87
  88       public:
  89         _Adaptor(_Engine& __g)
  90         : _M_g(__g) { }
  91
  92         _DInputType
  93         min() const
  94         { return _DInputType(0); }
  95
  96         _DInputType
  97         max() const
  98         { return _DInputType(1); }
  99
 100         /*
 101          * Converts a value generated by the adapted random number generator
 102          * into a value in the input domain for the dependent random number
 103          * distribution.
 104          */
 105         _DInputType
 106         operator()()
 107         {
 108           return std::generate_canonical<_DInputType,
 109                                     std::numeric_limits<_DInputType>::digits,
 110                                     _Engine>(_M_g);
 111         }
 112
 113       private:
 114         _Engine& _M_g;
 115       };
 116   } // namespace __detail
 117
 118   /**
 119    * @addtogroup random_generators Random Number Generators
 120    * @ingroup random
 121    *
 122    * These classes define objects which provide random or pseudorandom
 123    * numbers, either from a discrete or a continuous interval.  The
 124    * random number generator supplied as a part of this library are
 125    * all uniform random number generators which provide a sequence of
 126    * random number uniformly distributed over their range.
 127    *
 128    * A number generator is a function object with an operator() that
 129    * takes zero arguments and returns a number.
 130    *
 131    * A compliant random number generator must satisfy the following
 132    * requirements.  <table border=1 cellpadding=10 cellspacing=0>
 133    * <caption align=top>Random Number Generator Requirements</caption>
 134    * <tr><td>To be documented.</td></tr> </table>
 135    *
 136    * @{
 137    */
 138
 139   /**
 140    * @brief A model of a linear congruential random number generator.
 141    *
 142    * A random number generator that produces pseudorandom numbers via
 143    * linear function:
 144    * @f[
 145    *     x_{i+1}\leftarrow(ax_{i} + c) \bmod m
 146    * @f]
 147    *
 148    * The template parameter @p _UIntType must be an unsigned integral type
 149    * large enough to store values up to (__m-1). If the template parameter
 150    * @p __m is 0, the modulus @p __m used is
 151    * std::numeric_limits<_UIntType>::max() plus 1. Otherwise, the template
 152    * parameters @p __a and @p __c must be less than @p __m.
 153    *
 154    * The size of the state is @f$1@f$.
 155    */
 156   template<typename _UIntType, _UIntType __a, _UIntType __c, _UIntType __m>
 157     class linear_congruential_engine
 158     {
 159       static_assert(std::is_unsigned<_UIntType>::value, "template argument "
 160                     "substituting _UIntType not an unsigned integral type");
 161       static_assert(__m == 0u || (__a < __m && __c < __m),
 162                     "template argument substituting __m out of bounds");
 163
 164     public:
 165       /** The type of the generated random value. */
 166       typedef _UIntType result_type;
 167
 168       /** The multiplier. */
 169       static const result_type multiplier   = __a;
 170       /** An increment. */
 171       static const result_type increment    = __c;
 172       /** The modulus. */
 173       static const result_type modulus      = __m;
 174       static const result_type default_seed = 1u;
 175
 176       /**
 177        * @brief Constructs a %linear_congruential_engine random number
 178        *        generator engine with seed @p __s.  The default seed value
 179        *        is 1.
 180        *
 181        * @param __s The initial seed value.
 182        */
 183       explicit
 184       linear_congruential_engine(result_type __s = default_seed)
 185       { seed(__s); }
 186
 187       /**
 188        * @brief Constructs a %linear_congruential_engine random number
 189        *        generator engine seeded from the seed sequence @p __q.
 190        *
 191        * @param __q the seed sequence.
 192        */
 193       template<typename _Sseq, typename = typename
 194         std::enable_if<!std::is_same<_Sseq, linear_congruential_engine>::value>
 195                ::type>
 196         explicit
 197         linear_congruential_engine(_Sseq& __q)
 198         { seed(__q); }
 199
 200       /**
 201        * @brief Reseeds the %linear_congruential_engine random number generator
 202        *        engine sequence to the seed @p __s.
 203        *
 204        * @param __s The new seed.
 205        */
 206       void
 207       seed(result_type __s = default_seed);
 208
 209       /**
 210        * @brief Reseeds the %linear_congruential_engine random number generator
 211        *        engine
 212        * sequence using values from the seed sequence @p __q.
 213        *
 214        * @param __q the seed sequence.
 215        */
 216       template<typename _Sseq>
 217         typename std::enable_if<std::is_class<_Sseq>::value>::type
 218         seed(_Sseq& __q);
 219
 220       /**
 221        * @brief Gets the smallest possible value in the output range.
 222        *
 223        * The minimum depends on the @p __c parameter: if it is zero, the
 224        * minimum generated must be > 0, otherwise 0 is allowed.
 225        *
 226        * @todo This should be constexpr.
 227        */
 228       result_type
 229       min() const
 230       { return __c == 0u ? 1u : 0u; }
 231
 232       /**
 233        * @brief Gets the largest possible value in the output range.
 234        *
 235        * @todo This should be constexpr.
 236        */
 237       result_type
 238       max() const
 239       { return __m - 1u; }
 240
 241       /**
 242        * @brief Discard a sequence of random numbers.
 243        *
 244        * @todo Look for a faster way to do discard.
 245        */
 246       void
 247       discard(unsigned long long __z)
 248       {
 249         for (; __z != 0ULL; --__z)
 250           (*this)();
 251       }
 252
 253       /**
 254        * @brief Gets the next random number in the sequence.
 255        */
 256       result_type
 257       operator()()
 258       {
 259         _M_x = __detail::__mod<_UIntType, __m, __a, __c>(_M_x);
 260         return _M_x;
 261       }
 262
 263       /**
 264        * @brief Compares two linear congruential random number generator
 265        * objects of the same type for equality.
 266        *
 267        * @param __lhs A linear congruential random number generator object.
 268        * @param __rhs Another linear congruential random number generator
 269        *              object.
 270        *
 271        * @returns true if the two objects are equal, false otherwise.
 272        */
 273       friend bool
 274       operator==(const linear_congruential_engine& __lhs,
 275                  const linear_congruential_engine& __rhs)
 276       { return __lhs._M_x == __rhs._M_x; }
 277
 278       /**
 279        * @brief Writes the textual representation of the state x(i) of x to
 280        *        @p __os.
 281        *
 282        * @param __os  The output stream.
 283        * @param __lcr A % linear_congruential_engine random number generator.
 284        * @returns __os.
 285        */
 286       template<typename _UIntType1, _UIntType1 __a1, _UIntType1 __c1,
 287                _UIntType1 __m1, typename _CharT, typename _Traits>
 288         friend std::basic_ostream<_CharT, _Traits>&
 289         operator<<(std::basic_ostream<_CharT, _Traits>&,
 290                    const std::linear_congruential_engine<_UIntType1,
 291                    __a1, __c1, __m1>&);
 292
 293       /**
 294        * @brief Sets the state of the engine by reading its textual
 295        *        representation from @p __is.
 296        *
 297        * The textual representation must have been previously written using
 298        * an output stream whose imbued locale and whose type's template
 299        * specialization arguments _CharT and _Traits were the same as those
 300        * of @p __is.
 301        *
 302        * @param __is  The input stream.
 303        * @param __lcr A % linear_congruential_engine random number generator.
 304        * @returns __is.
 305        */
 306       template<typename _UIntType1, _UIntType1 __a1, _UIntType1 __c1,
 307                _UIntType1 __m1, typename _CharT, typename _Traits>
 308         friend std::basic_istream<_CharT, _Traits>&
 309         operator>>(std::basic_istream<_CharT, _Traits>&,
 310                    std::linear_congruential_engine<_UIntType1, __a1,
 311                    __c1, __m1>&);
 312
 313     private:
 314       _UIntType _M_x;
 315     };
 316
 317
 318   /**
 319    * A generalized feedback shift register discrete random number generator.
 320    *
 321    * This algorithm avoids multiplication and division and is designed to be
 322    * friendly to a pipelined architecture.  If the parameters are chosen
 323    * correctly, this generator will produce numbers with a very long period and
 324    * fairly good apparent entropy, although still not cryptographically strong.
 325    *
 326    * The best way to use this generator is with the predefined mt19937 class.
 327    *
 328    * This algorithm was originally invented by Makoto Matsumoto and
 329    * Takuji Nishimura.
 330    *
 331    * @var word_size   The number of bits in each element of the state vector.
 332    * @var state_size  The degree of recursion.
 333    * @var shift_size  The period parameter.
 334    * @var mask_bits   The separation point bit index.
 335    * @var parameter_a The last row of the twist matrix.
 336    * @var output_u    The first right-shift tempering matrix parameter.
 337    * @var output_s    The first left-shift tempering matrix parameter.
 338    * @var output_b    The first left-shift tempering matrix mask.
 339    * @var output_t    The second left-shift tempering matrix parameter.
 340    * @var output_c    The second left-shift tempering matrix mask.
 341    * @var output_l    The second right-shift tempering matrix parameter.
 342    */
 343   template<typename _UIntType, size_t __w,
 344            size_t __n, size_t __m, size_t __r,
 345            _UIntType __a, size_t __u, _UIntType __d, size_t __s,
 346            _UIntType __b, size_t __t,
 347            _UIntType __c, size_t __l, _UIntType __f>
 348     class mersenne_twister_engine
 349     {
 350       static_assert(std::is_unsigned<_UIntType>::value, "template argument "
 351                     "substituting _UIntType not an unsigned integral type");
 352       static_assert(1u <= __m && __m <= __n,
 353                     "template argument substituting __m out of bounds");
 354       static_assert(__r <= __w, "template argument substituting "
 355                     "__r out of bound");
 356       static_assert(__u <= __w, "template argument substituting "
 357                     "__u out of bound");
 358       static_assert(__s <= __w, "template argument substituting "
 359                     "__s out of bound");
 360       static_assert(__t <= __w, "template argument substituting "
 361                     "__t out of bound");
 362       static_assert(__l <= __w, "template argument substituting "
 363                     "__l out of bound");
 364       static_assert(__w <= std::numeric_limits<_UIntType>::digits,
 365                     "template argument substituting __w out of bound");
 366       static_assert(__a <= (__detail::_Shift<_UIntType, __w>::__value - 1),
 367                     "template argument substituting __a out of bound");
 368       static_assert(__b <= (__detail::_Shift<_UIntType, __w>::__value - 1),
 369                     "template argument substituting __b out of bound");
 370       static_assert(__c <= (__detail::_Shift<_UIntType, __w>::__value - 1),
 371                     "template argument substituting __c out of bound");
 372       static_assert(__d <= (__detail::_Shift<_UIntType, __w>::__value - 1),
 373                     "template argument substituting __d out of bound");
 374       static_assert(__f <= (__detail::_Shift<_UIntType, __w>::__value - 1),
 375                     "template argument substituting __f out of bound");
 376
 377     public:
 378       /** The type of the generated random value. */
 379       typedef _UIntType result_type;
 380
 381       // parameter values
 382       static const size_t      word_size                 = __w;
 383       static const size_t      state_size                = __n;
 384       static const size_t      shift_size                = __m;
 385       static const size_t      mask_bits                 = __r;
 386       static const result_type xor_mask                  = __a;
 387       static const size_t      tempering_u               = __u;
 388       static const result_type tempering_d               = __d;
 389       static const size_t      tempering_s               = __s;
 390       static const result_type tempering_b               = __b;
 391       static const size_t      tempering_t               = __t;
 392       static const result_type tempering_c               = __c;
 393       static const size_t      tempering_l               = __l;
 394       static const result_type initialization_multiplier = __f;
 395       static const result_type default_seed = 5489u;
 396
 397       // constructors and member function
 398       explicit
 399       mersenne_twister_engine(result_type __sd = default_seed)
 400       { seed(__sd); }
 401
 402       /**
 403        * @brief Constructs a %mersenne_twister_engine random number generator
 404        *        engine seeded from the seed sequence @p __q.
 405        *
 406        * @param __q the seed sequence.
 407        */
 408       template<typename _Sseq, typename = typename
 409         std::enable_if<!std::is_same<_Sseq, mersenne_twister_engine>::value>
 410                ::type>
 411         explicit
 412         mersenne_twister_engine(_Sseq& __q)
 413         { seed(__q); }
 414
 415       void
 416       seed(result_type __sd = default_seed);
 417
 418       template<typename _Sseq>
 419         typename std::enable_if<std::is_class<_Sseq>::value>::type
 420         seed(_Sseq& __q);
 421
 422       /**
 423        * @brief Gets the smallest possible value in the output range.
 424        *
 425        * @todo This should be constexpr.
 426        */
 427       result_type
 428       min() const
 429       { return 0; };
 430
 431       /**
 432        * @brief Gets the largest possible value in the output range.
 433        *
 434        * @todo This should be constexpr.
 435        */
 436       result_type
 437       max() const
 438       { return __detail::_Shift<_UIntType, __w>::__value - 1; }
 439
 440       /**
 441        * @brief Discard a sequence of random numbers.
 442        *
 443        * @todo Look for a faster way to do discard.
 444        */
 445       void
 446       discard(unsigned long long __z)
 447       {
 448         for (; __z != 0ULL; --__z)
 449           (*this)();
 450       }
 451
 452       result_type
 453       operator()();
 454
 455       /**
 456        * @brief Compares two % mersenne_twister_engine random number generator
 457        *        objects of the same type for equality.
 458        *
 459        * @param __lhs A % mersenne_twister_engine random number generator
 460        *              object.
 461        * @param __rhs Another % mersenne_twister_engine random number
 462        *              generator object.
 463        *
 464        * @returns true if the two objects are equal, false otherwise.
 465        */
 466       friend bool
 467       operator==(const mersenne_twister_engine& __lhs,
 468                  const mersenne_twister_engine& __rhs)
 469       { return std::equal(__lhs._M_x, __lhs._M_x + state_size, __rhs._M_x); }
 470
 471       /**
 472        * @brief Inserts the current state of a % mersenne_twister_engine
 473        *        random number generator engine @p __x into the output stream
 474        *        @p __os.
 475        *
 476        * @param __os An output stream.
 477        * @param __x  A % mersenne_twister_engine random number generator
 478        *             engine.
 479        *
 480        * @returns The output stream with the state of @p __x inserted or in
 481        * an error state.
 482        */
 483       template<typename _UIntType1,
 484                size_t __w1, size_t __n1,
 485                size_t __m1, size_t __r1,
 486                _UIntType1 __a1, size_t __u1,
 487                _UIntType1 __d1, size_t __s1,
 488                _UIntType1 __b1, size_t __t1,
 489                _UIntType1 __c1, size_t __l1, _UIntType1 __f1,
 490                typename _CharT, typename _Traits>
 491         friend std::basic_ostream<_CharT, _Traits>&
 492         operator<<(std::basic_ostream<_CharT, _Traits>&,
 493                    const std::mersenne_twister_engine<_UIntType1, __w1, __n1,
 494                    __m1, __r1, __a1, __u1, __d1, __s1, __b1, __t1, __c1,
 495                    __l1, __f1>&);
 496
 497       /**
 498        * @brief Extracts the current state of a % mersenne_twister_engine
 499        *        random number generator engine @p __x from the input stream
 500        *        @p __is.
 501        *
 502        * @param __is An input stream.
 503        * @param __x  A % mersenne_twister_engine random number generator
 504        *             engine.
 505        *
 506        * @returns The input stream with the state of @p __x extracted or in
 507        * an error state.
 508        */
 509       template<typename _UIntType1,
 510                size_t __w1, size_t __n1,
 511                size_t __m1, size_t __r1,
 512                _UIntType1 __a1, size_t __u1,
 513                _UIntType1 __d1, size_t __s1,
 514                _UIntType1 __b1, size_t __t1,
 515                _UIntType1 __c1, size_t __l1, _UIntType1 __f1,
 516                typename _CharT, typename _Traits>
 517         friend std::basic_istream<_CharT, _Traits>&
 518         operator>>(std::basic_istream<_CharT, _Traits>&,
 519                    std::mersenne_twister_engine<_UIntType1, __w1, __n1, __m1,
 520                    __r1, __a1, __u1, __d1, __s1, __b1, __t1, __c1,
 521                    __l1, __f1>&);
 522
 523     private:
 524       _UIntType _M_x[state_size];
 525       size_t    _M_p;
 526     };
 527
 528   /**
 529    * @brief The Marsaglia-Zaman generator.
 530    *
 531    * This is a model of a Generalized Fibonacci discrete random number
 532    * generator, sometimes referred to as the SWC generator.
 533    *
 534    * A discrete random number generator that produces pseudorandom
 535    * numbers using:
 536    * @f[
 537    *     x_{i}\leftarrow(x_{i - s} - x_{i - r} - carry_{i-1}) \bmod m
 538    * @f]
 539    *
 540    * The size of the state is @f$r@f$
 541    * and the maximum period of the generator is @f$(m^r - m^s - 1)@f$.
 542    *
 543    * @var _M_x     The state of the generator.  This is a ring buffer.
 544    * @var _M_carry The carry.
 545    * @var _M_p     Current index of x(i - r).
 546    */
 547   template<typename _UIntType, size_t __w, size_t __s, size_t __r>
 548     class subtract_with_carry_engine
 549     {
 550       static_assert(std::is_unsigned<_UIntType>::value, "template argument "
 551                     "substituting _UIntType not an unsigned integral type");
 552       static_assert(0u < __s && __s < __r,
 553                     "template argument substituting __s out of bounds");
 554       static_assert(0u < __w && __w <= std::numeric_limits<_UIntType>::digits,
 555                     "template argument substituting __w out of bounds");
 556
 557     public:
 558       /** The type of the generated random value. */
 559       typedef _UIntType result_type;
 560
 561       // parameter values
 562       static const size_t      word_size    = __w;
 563       static const size_t      short_lag    = __s;
 564       static const size_t      long_lag     = __r;
 565       static const result_type default_seed = 19780503u;
 566
 567       /**
 568        * @brief Constructs an explicitly seeded % subtract_with_carry_engine
 569        *        random number generator.
 570        */
 571       explicit
 572       subtract_with_carry_engine(result_type __sd = default_seed)
 573       { seed(__sd); }
 574
 575       /**
 576        * @brief Constructs a %subtract_with_carry_engine random number engine
 577        *        seeded from the seed sequence @p __q.
 578        *
 579        * @param __q the seed sequence.
 580        */
 581       template<typename _Sseq, typename = typename
 582         std::enable_if<!std::is_same<_Sseq, subtract_with_carry_engine>::value>
 583                ::type>
 584         explicit
 585         subtract_with_carry_engine(_Sseq& __q)
 586         { seed(__q); }
 587
 588       /**
 589        * @brief Seeds the initial state @f$x_0@f$ of the random number
 590        *        generator.
 591        *
 592        * N1688[4.19] modifies this as follows.  If @p __value == 0,
 593        * sets value to 19780503.  In any case, with a linear
 594        * congruential generator lcg(i) having parameters @f$ m_{lcg} =
 595        * 2147483563, a_{lcg} = 40014, c_{lcg} = 0, and lcg(0) = value
 596        * @f$, sets @f$ x_{-r} \dots x_{-1} @f$ to @f$ lcg(1) \bmod m
 597        * \dots lcg(r) \bmod m @f$ respectively.  If @f$ x_{-1} = 0 @f$
 598        * set carry to 1, otherwise sets carry to 0.
 599        */
 600       void
 601       seed(result_type __sd = default_seed);
 602
 603       /**
 604        * @brief Seeds the initial state @f$x_0@f$ of the
 605        * % subtract_with_carry_engine random number generator.
 606        */
 607       template<typename _Sseq>
 608         typename std::enable_if<std::is_class<_Sseq>::value>::type
 609         seed(_Sseq& __q);
 610
 611       /**
 612        * @brief Gets the inclusive minimum value of the range of random
 613        * integers returned by this generator.
 614        *
 615        * @todo This should be constexpr.
 616        */
 617       result_type
 618       min() const
 619       { return 0; }
 620
 621       /**
 622        * @brief Gets the inclusive maximum value of the range of random
 623        * integers returned by this generator.
 624        *
 625        * @todo This should be constexpr.
 626        */
 627       result_type
 628       max() const
 629       { return __detail::_Shift<_UIntType, __w>::__value - 1; }
 630
 631       /**
 632        * @brief Discard a sequence of random numbers.
 633        *
 634        * @todo Look for a faster way to do discard.
 635        */
 636       void
 637       discard(unsigned long long __z)
 638       {
 639         for (; __z != 0ULL; --__z)
 640           (*this)();
 641       }
 642
 643       /**
 644        * @brief Gets the next random number in the sequence.
 645        */
 646       result_type
 647       operator()();
 648
 649       /**
 650        * @brief Compares two % subtract_with_carry_engine random number
 651        *        generator objects of the same type for equality.
 652        *
 653        * @param __lhs A % subtract_with_carry_engine random number generator
 654        *              object.
 655        * @param __rhs Another % subtract_with_carry_engine random number
 656        *              generator object.
 657        *
 658        * @returns true if the two objects are equal, false otherwise.
 659        */
 660       friend bool
 661       operator==(const subtract_with_carry_engine& __lhs,
 662                  const subtract_with_carry_engine& __rhs)
 663       { return std::equal(__lhs._M_x, __lhs._M_x + long_lag, __rhs._M_x); }
 664
 665       /**
 666        * @brief Inserts the current state of a % subtract_with_carry_engine
 667        *        random number generator engine @p __x into the output stream
 668        *        @p __os.
 669        *
 670        * @param __os An output stream.
 671        * @param __x  A % subtract_with_carry_engine random number generator
 672        *             engine.
 673        *
 674        * @returns The output stream with the state of @p __x inserted or in
 675        * an error state.
 676        */
 677       template<typename _UIntType1, size_t __w1, size_t __s1, size_t __r1,
 678                typename _CharT, typename _Traits>
 679         friend std::basic_ostream<_CharT, _Traits>&
 680         operator<<(std::basic_ostream<_CharT, _Traits>&,
 681                    const std::subtract_with_carry_engine<_UIntType1, __w1,
 682                    __s1, __r1>&);
 683
 684       /**
 685        * @brief Extracts the current state of a % subtract_with_carry_engine
 686        *        random number generator engine @p __x from the input stream
 687        *        @p __is.
 688        *
 689        * @param __is An input stream.
 690        * @param __x  A % subtract_with_carry_engine random number generator
 691        *             engine.
 692        *
 693        * @returns The input stream with the state of @p __x extracted or in
 694        * an error state.
 695        */
 696       template<typename _UIntType1, size_t __w1, size_t __s1, size_t __r1,
 697                typename _CharT, typename _Traits>
 698         friend std::basic_istream<_CharT, _Traits>&
 699         operator>>(std::basic_istream<_CharT, _Traits>&,
 700                    std::subtract_with_carry_engine<_UIntType1, __w1,
 701                    __s1, __r1>&);
 702
 703     private:
 704       _UIntType  _M_x[long_lag];
 705       _UIntType  _M_carry;
 706       size_t     _M_p;
 707     };
 708
 709   /**
 710    * Produces random numbers from some base engine by discarding blocks of
 711    * data.
 712    *
 713    * 0 <= @p __r <= @p __p
 714    */
 715   template<typename _RandomNumberEngine, size_t __p, size_t __r>
 716     class discard_block_engine
 717     {
 718       static_assert(1 <= __r && __r <= __p,
 719                     "template argument substituting __r out of bounds");
 720
 721     public:
 722       /** The type of the generated random value. */
 723       typedef typename _RandomNumberEngine::result_type result_type;
 724
 725       // parameter values
 726       static const size_t block_size = __p;
 727       static const size_t used_block = __r;
 728
 729       /**
 730        * @brief Constructs a default %discard_block_engine engine.
 731        *
 732        * The underlying engine is default constructed as well.
 733        */
 734       discard_block_engine()
 735       : _M_b(), _M_n(0) { }
 736
 737       /**
 738        * @brief Copy constructs a %discard_block_engine engine.
 739        *
 740        * Copies an existing base class random number generator.
 741        * @param rng An existing (base class) engine object.
 742        */
 743       explicit
 744       discard_block_engine(const _RandomNumberEngine& __rne)
 745       : _M_b(__rne), _M_n(0) { }
 746
 747       /**
 748        * @brief Move constructs a %discard_block_engine engine.
 749        *
 750        * Copies an existing base class random number generator.
 751        * @param rng An existing (base class) engine object.
 752        */
 753       explicit
 754       discard_block_engine(_RandomNumberEngine&& __rne)
 755       : _M_b(std::move(__rne)), _M_n(0) { }
 756
 757       /**
 758        * @brief Seed constructs a %discard_block_engine engine.
 759        *
 760        * Constructs the underlying generator engine seeded with @p __s.
 761        * @param __s A seed value for the base class engine.
 762        */
 763       explicit
 764       discard_block_engine(result_type __s)
 765       : _M_b(__s), _M_n(0) { }
 766
 767       /**
 768        * @brief Generator construct a %discard_block_engine engine.
 769        *
 770        * @param __q A seed sequence.
 771        */
 772       template<typename _Sseq, typename = typename
 773         std::enable_if<!std::is_same<_Sseq, discard_block_engine>::value
 774                        && !std::is_same<_Sseq, _RandomNumberEngine>::value>
 775                ::type>
 776         explicit
 777         discard_block_engine(_Sseq& __q)
 778         : _M_b(__q), _M_n(0)
 779         { }
 780
 781       /**
 782        * @brief Reseeds the %discard_block_engine object with the default
 783        *        seed for the underlying base class generator engine.
 784        */
 785       void
 786       seed()
 787       {
 788         _M_b.seed();
 789         _M_n = 0;
 790       }
 791
 792       /**
 793        * @brief Reseeds the %discard_block_engine object with the default
 794        *        seed for the underlying base class generator engine.
 795        */
 796       void
 797       seed(result_type __s)
 798       {
 799         _M_b.seed(__s);
 800         _M_n = 0;
 801       }
 802
 803       /**
 804        * @brief Reseeds the %discard_block_engine object with the given seed
 805        *        sequence.
 806        * @param __q A seed generator function.
 807        */
 808       template<typename _Sseq>
 809         void
 810         seed(_Sseq& __q)
 811         {
 812           _M_b.seed(__q);
 813           _M_n = 0;
 814         }
 815
 816       /**
 817        * @brief Gets a const reference to the underlying generator engine
 818        *        object.
 819        */
 820       const _RandomNumberEngine&
 821       base() const
 822       { return _M_b; }
 823
 824       /**
 825        * @brief Gets the minimum value in the generated random number range.
 826        *
 827        * @todo This should be constexpr.
 828        */
 829       result_type
 830       min() const
 831       { return _M_b.min(); }
 832
 833       /**
 834        * @brief Gets the maximum value in the generated random number range.
 835        *
 836        * @todo This should be constexpr.
 837        */
 838       result_type
 839       max() const
 840       { return _M_b.max(); }
 841
 842       /**
 843        * @brief Discard a sequence of random numbers.
 844        *
 845        * @todo Look for a faster way to do discard.
 846        */
 847       void
 848       discard(unsigned long long __z)
 849       {
 850         for (; __z != 0ULL; --__z)
 851           (*this)();
 852       }
 853
 854       /**
 855        * @brief Gets the next value in the generated random number sequence.
 856        */
 857       result_type
 858       operator()();
 859
 860       /**
 861        * @brief Compares two %discard_block_engine random number generator
 862        *        objects of the same type for equality.
 863        *
 864        * @param __lhs A %discard_block_engine random number generator object.
 865        * @param __rhs Another %discard_block_engine random number generator
 866        *              object.
 867        *
 868        * @returns true if the two objects are equal, false otherwise.
 869        */
 870       friend bool
 871       operator==(const discard_block_engine& __lhs,
 872                  const discard_block_engine& __rhs)
 873       { return (__lhs._M_b == __rhs._M_b) && (__lhs._M_n == __rhs._M_n); }
 874
 875       /**
 876        * @brief Inserts the current state of a %discard_block_engine random
 877        *        number generator engine @p __x into the output stream
 878        *        @p __os.
 879        *
 880        * @param __os An output stream.
 881        * @param __x  A %discard_block_engine random number generator engine.
 882        *
 883        * @returns The output stream with the state of @p __x inserted or in
 884        * an error state.
 885        */
 886       template<typename _RandomNumberEngine1, size_t __p1, size_t __r1,
 887                typename _CharT, typename _Traits>
 888         friend std::basic_ostream<_CharT, _Traits>&
 889         operator<<(std::basic_ostream<_CharT, _Traits>&,
 890                    const std::discard_block_engine<_RandomNumberEngine1,
 891                    __p1, __r1>&);
 892
 893       /**
 894        * @brief Extracts the current state of a % subtract_with_carry_engine
 895        *        random number generator engine @p __x from the input stream
 896        *        @p __is.
 897        *
 898        * @param __is An input stream.
 899        * @param __x  A %discard_block_engine random number generator engine.
 900        *
 901        * @returns The input stream with the state of @p __x extracted or in
 902        * an error state.
 903        */
 904       template<typename _RandomNumberEngine1, size_t __p1, size_t __r1,
 905                typename _CharT, typename _Traits>
 906         friend std::basic_istream<_CharT, _Traits>&
 907         operator>>(std::basic_istream<_CharT, _Traits>&,
 908                    std::discard_block_engine<_RandomNumberEngine1,
 909                    __p1, __r1>&);
 910
 911     private:
 912       _RandomNumberEngine _M_b;
 913       size_t _M_n;
 914     };
 915
 916   /**
 917    * Produces random numbers by combining random numbers from some base
 918    * engine to produce random numbers with a specifies number of bits @p __w.
 919    */
 920   template<typename _RandomNumberEngine, size_t __w, typename _UIntType>
 921     class independent_bits_engine
 922     {
 923       static_assert(std::is_unsigned<_UIntType>::value, "template argument "
 924                     "substituting _UIntType not an unsigned integral type");
 925       static_assert(0u < __w && __w <= std::numeric_limits<_UIntType>::digits,
 926                     "template argument substituting __w out of bounds");
 927
 928     public:
 929       /** The type of the generated random value. */
 930       typedef _UIntType result_type;
 931
 932       /**
 933        * @brief Constructs a default %independent_bits_engine engine.
 934        *
 935        * The underlying engine is default constructed as well.
 936        */
 937       independent_bits_engine()
 938       : _M_b() { }
 939
 940       /**
 941        * @brief Copy constructs a %independent_bits_engine engine.
 942        *
 943        * Copies an existing base class random number generator.
 944        * @param rng An existing (base class) engine object.
 945        */
 946       explicit
 947       independent_bits_engine(const _RandomNumberEngine& __rne)
 948       : _M_b(__rne) { }
 949
 950       /**
 951        * @brief Move constructs a %independent_bits_engine engine.
 952        *
 953        * Copies an existing base class random number generator.
 954        * @param rng An existing (base class) engine object.
 955        */
 956       explicit
 957       independent_bits_engine(_RandomNumberEngine&& __rne)
 958       : _M_b(std::move(__rne)) { }
 959
 960       /**
 961        * @brief Seed constructs a %independent_bits_engine engine.
 962        *
 963        * Constructs the underlying generator engine seeded with @p __s.
 964        * @param __s A seed value for the base class engine.
 965        */
 966       explicit
 967       independent_bits_engine(result_type __s)
 968       : _M_b(__s) { }
 969
 970       /**
 971        * @brief Generator construct a %independent_bits_engine engine.
 972        *
 973        * @param __q A seed sequence.
 974        */
 975       template<typename _Sseq, typename = typename
 976         std::enable_if<!std::is_same<_Sseq, independent_bits_engine>::value
 977                        && !std::is_same<_Sseq, _RandomNumberEngine>::value>
 978                ::type>
 979         explicit
 980         independent_bits_engine(_Sseq& __q)
 981         : _M_b(__q)
 982         { }
 983
 984       /**
 985        * @brief Reseeds the %independent_bits_engine object with the default
 986        *        seed for the underlying base class generator engine.
 987        */
 988       void
 989       seed()
 990       { _M_b.seed(); }
 991
 992       /**
 993        * @brief Reseeds the %independent_bits_engine object with the default
 994        *        seed for the underlying base class generator engine.
 995        */
 996       void
 997       seed(result_type __s)
 998       { _M_b.seed(__s); }
 999
1000       /**
1001        * @brief Reseeds the %independent_bits_engine object with the given
1002        *        seed sequence.
1003        * @param __q A seed generator function.
1004        */
1005       template<typename _Sseq>
1006         void
1007         seed(_Sseq& __q)
1008         { _M_b.seed(__q); }
1009
1010       /**
1011        * @brief Gets a const reference to the underlying generator engine
1012        *        object.
1013        */
1014       const _RandomNumberEngine&
1015       base() const
1016       { return _M_b; }
1017
1018       /**
1019        * @brief Gets the minimum value in the generated random number range.
1020        *
1021        * @todo This should be constexpr.
1022        */
1023       result_type
1024       min() const
1025       { return 0U; }
1026
1027       /**
1028        * @brief Gets the maximum value in the generated random number range.
1029        *
1030        * @todo This should be constexpr.
1031        */
1032       result_type
1033       max() const
1034       { return __detail::_Shift<_UIntType, __w>::__value - 1; }
1035
1036       /**
1037        * @brief Discard a sequence of random numbers.
1038        *
1039        * @todo Look for a faster way to do discard.
1040        */
1041       void
1042       discard(unsigned long long __z)
1043       {
1044         for (; __z != 0ULL; --__z)
1045           (*this)();
1046       }
1047
1048       /**
1049        * @brief Gets the next value in the generated random number sequence.
1050        */
1051       result_type
1052       operator()();
1053
1054       /**
1055        * @brief Compares two %independent_bits_engine random number generator
1056        * objects of the same type for equality.
1057        *
1058        * @param __lhs A %independent_bits_engine random number generator
1059        *              object.
1060        * @param __rhs Another %independent_bits_engine random number generator
1061        *              object.
1062        *
1063        * @returns true if the two objects are equal, false otherwise.
1064        */
1065       friend bool
1066       operator==(const independent_bits_engine& __lhs,
1067                  const independent_bits_engine& __rhs)
1068       { return __lhs._M_b == __rhs._M_b; }
1069
1070       /**
1071        * @brief Extracts the current state of a % subtract_with_carry_engine
1072        *        random number generator engine @p __x from the input stream
1073        *        @p __is.
1074        *
1075        * @param __is An input stream.
1076        * @param __x  A %independent_bits_engine random number generator
1077        *             engine.
1078        *
1079        * @returns The input stream with the state of @p __x extracted or in
1080        *          an error state.
1081        */
1082       template<typename _CharT, typename _Traits>
1083         friend std::basic_istream<_CharT, _Traits>&
1084         operator>>(std::basic_istream<_CharT, _Traits>& __is,
1085                    std::independent_bits_engine<_RandomNumberEngine,
1086                    __w, _UIntType>& __x)
1087         {
1088           __is >> __x._M_b;
1089           return __is;
1090         }
1091
1092     private:
1093       _RandomNumberEngine _M_b;
1094     };
1095
1096   /**
1097    * @brief Inserts the current state of a %independent_bits_engine random
1098    *        number generator engine @p __x into the output stream @p __os.
1099    *
1100    * @param __os An output stream.
1101    * @param __x  A %independent_bits_engine random number generator engine.
1102    *
1103    * @returns The output stream with the state of @p __x inserted or in
1104    *          an error state.
1105    */
1106   template<typename _RandomNumberEngine, size_t __w, typename _UIntType,
1107            typename _CharT, typename _Traits>
1108     std::basic_ostream<_CharT, _Traits>&
1109     operator<<(std::basic_ostream<_CharT, _Traits>& __os,
1110                const std::independent_bits_engine<_RandomNumberEngine,
1111                __w, _UIntType>& __x)
1112     {
1113       __os << __x.base();
1114       return __os;
1115     }
1116
1117   /**
1118    * @brief Produces random numbers by combining random numbers from some
1119    * base engine to produce random numbers with a specifies number of bits
1120    * @p __w.
1121    */
1122   template<typename _RandomNumberEngine, size_t __k>
1123     class shuffle_order_engine
1124     {
1125       static_assert(1u <= __k, "template argument substituting "
1126                     "__k out of bound");
1127
1128     public:
1129       /** The type of the generated random value. */
1130       typedef typename _RandomNumberEngine::result_type result_type;
1131
1132       static const size_t table_size = __k;
1133
1134       /**
1135        * @brief Constructs a default %shuffle_order_engine engine.
1136        *
1137        * The underlying engine is default constructed as well.
1138        */
1139       shuffle_order_engine()
1140       : _M_b()
1141       { _M_initialize(); }
1142
1143       /**
1144        * @brief Copy constructs a %shuffle_order_engine engine.
1145        *
1146        * Copies an existing base class random number generator.
1147        * @param rng An existing (base class) engine object.
1148        */
1149       explicit
1150       shuffle_order_engine(const _RandomNumberEngine& __rne)
1151       : _M_b(__rne)
1152       { _M_initialize(); }
1153
1154       /**
1155        * @brief Move constructs a %shuffle_order_engine engine.
1156        *
1157        * Copies an existing base class random number generator.
1158        * @param rng An existing (base class) engine object.
1159        */
1160       explicit
1161       shuffle_order_engine(_RandomNumberEngine&& __rne)
1162       : _M_b(std::move(__rne))
1163       { _M_initialize(); }
1164
1165       /**
1166        * @brief Seed constructs a %shuffle_order_engine engine.
1167        *
1168        * Constructs the underlying generator engine seeded with @p __s.
1169        * @param __s A seed value for the base class engine.
1170        */
1171       explicit
1172       shuffle_order_engine(result_type __s)
1173       : _M_b(__s)
1174       { _M_initialize(); }
1175
1176       /**
1177        * @brief Generator construct a %shuffle_order_engine engine.
1178        *
1179        * @param __q A seed sequence.
1180        */
1181       template<typename _Sseq, typename = typename
1182         std::enable_if<!std::is_same<_Sseq, shuffle_order_engine>::value
1183                        && !std::is_same<_Sseq, _RandomNumberEngine>::value>
1184                ::type>
1185         explicit
1186         shuffle_order_engine(_Sseq& __q)
1187         : _M_b(__q)
1188         { _M_initialize(); }
1189
1190       /**
1191        * @brief Reseeds the %shuffle_order_engine object with the default seed
1192                 for the underlying base class generator engine.
1193        */
1194       void
1195       seed()
1196       {
1197         _M_b.seed();
1198         _M_initialize();
1199       }
1200
1201       /**
1202        * @brief Reseeds the %shuffle_order_engine object with the default seed
1203        *        for the underlying base class generator engine.
1204        */
1205       void
1206       seed(result_type __s)
1207       {
1208         _M_b.seed(__s);
1209         _M_initialize();
1210       }
1211
1212       /**
1213        * @brief Reseeds the %shuffle_order_engine object with the given seed
1214        *        sequence.
1215        * @param __q A seed generator function.
1216        */
1217       template<typename _Sseq>
1218         void
1219         seed(_Sseq& __q)
1220         {
1221           _M_b.seed(__q);
1222           _M_initialize();
1223         }
1224
1225       /**
1226        * Gets a const reference to the underlying generator engine object.
1227        */
1228       const _RandomNumberEngine&
1229       base() const
1230       { return _M_b; }
1231
1232       /**
1233        * Gets the minimum value in the generated random number range.
1234        *
1235        * @todo This should be constexpr.
1236        */
1237       result_type
1238       min() const
1239       { return _M_b.min(); }
1240
1241       /**
1242        * Gets the maximum value in the generated random number range.
1243        *
1244        * @todo This should be constexpr.
1245        */
1246       result_type
1247       max() const
1248       { return _M_b.max(); }
1249
1250       /**
1251        * Discard a sequence of random numbers.
1252        *
1253        * @todo Look for a faster way to do discard.
1254        */
1255       void
1256       discard(unsigned long long __z)
1257       {
1258         for (; __z != 0ULL; --__z)
1259           (*this)();
1260       }
1261
1262       /**
1263        * Gets the next value in the generated random number sequence.
1264        */
1265       result_type
1266       operator()();
1267
1268       /**
1269        * Compares two %shuffle_order_engine random number generator objects
1270        * of the same type for equality.
1271        *
1272        * @param __lhs A %shuffle_order_engine random number generator object.
1273        * @param __rhs Another %shuffle_order_engine random number generator
1274        *              object.
1275        *
1276        * @returns true if the two objects are equal, false otherwise.
1277        */
1278       friend bool
1279       operator==(const shuffle_order_engine& __lhs,
1280                  const shuffle_order_engine& __rhs)
1281       { return __lhs._M_b == __rhs._M_b; }
1282
1283       /**
1284        * @brief Inserts the current state of a %shuffle_order_engine random
1285        *        number generator engine @p __x into the output stream
1286         @p __os.
1287        *
1288        * @param __os An output stream.
1289        * @param __x  A %shuffle_order_engine random number generator engine.
1290        *
1291        * @returns The output stream with the state of @p __x inserted or in
1292        * an error state.
1293        */
1294       template<typename _RandomNumberEngine1, size_t __k1,
1295                typename _CharT, typename _Traits>
1296         friend std::basic_ostream<_CharT, _Traits>&
1297         operator<<(std::basic_ostream<_CharT, _Traits>&,
1298                    const std::shuffle_order_engine<_RandomNumberEngine1,
1299                    __k1>&);
1300
1301       /**
1302        * @brief Extracts the current state of a % subtract_with_carry_engine
1303        *        random number generator engine @p __x from the input stream
1304        *        @p __is.
1305        *
1306        * @param __is An input stream.
1307        * @param __x  A %shuffle_order_engine random number generator engine.
1308        *
1309        * @returns The input stream with the state of @p __x extracted or in
1310        * an error state.
1311        */
1312       template<typename _RandomNumberEngine1, size_t __k1,
1313                typename _CharT, typename _Traits>
1314         friend std::basic_istream<_CharT, _Traits>&
1315         operator>>(std::basic_istream<_CharT, _Traits>&,
1316                    std::shuffle_order_engine<_RandomNumberEngine1, __k1>&);
1317
1318     private:
1319       void _M_initialize()
1320       {
1321         for (size_t __i = 0; __i < __k; ++__i)
1322           _M_v[__i] = _M_b();
1323         _M_y = _M_b();
1324       }
1325
1326       _RandomNumberEngine _M_b;
1327       result_type _M_v[__k];
1328       result_type _M_y;
1329     };
1330
1331   /**
1332    * The classic Minimum Standard rand0 of Lewis, Goodman, and Miller.
1333    */
1334   typedef linear_congruential_engine<uint_fast32_t, 16807UL, 0UL, 2147483647UL>
1335   minstd_rand0;
1336
1337   /**
1338    * An alternative LCR (Lehmer Generator function).
1339    */
1340   typedef linear_congruential_engine<uint_fast32_t, 48271UL, 0UL, 2147483647UL>
1341   minstd_rand;
1342
1343   /**
1344    * The classic Mersenne Twister.
1345    *
1346    * Reference:
1347    * M. Matsumoto and T. Nishimura, Mersenne Twister: A 623-Dimensionally
1348    * Equidistributed Uniform Pseudo-Random Number Generator, ACM Transactions
1349    * on Modeling and Computer Simulation, Vol. 8, No. 1, January 1998, pp 3-30.
1350    */
1351   typedef mersenne_twister_engine<
1352     uint_fast32_t,
1353     32, 624, 397, 31,
1354     0x9908b0dfUL, 11,
1355     0xffffffffUL, 7,
1356     0x9d2c5680UL, 15,
1357     0xefc60000UL, 18, 1812433253UL> mt19937;
1358
1359   /**
1360    * An alternative Mersenne Twister.
1361    */
1362   typedef mersenne_twister_engine<
1363     uint_fast64_t,
1364     64, 312, 156, 31,
1365     0xb5026f5aa96619e9ULL, 29,
1366     0x5555555555555555ULL, 17,
1367     0x71d67fffeda60000ULL, 37,
1368     0xfff7eee000000000ULL, 43,
1369     6364136223846793005ULL> mt19937_64;
1370
1371   typedef subtract_with_carry_engine<uint_fast32_t, 24, 10, 24>
1372     ranlux24_base;
1373
1374   typedef subtract_with_carry_engine<uint_fast64_t, 48, 5, 12>
1375     ranlux48_base;
1376
1377   typedef discard_block_engine<ranlux24_base, 223, 23> ranlux24;
1378
1379   typedef discard_block_engine<ranlux48_base, 389, 11> ranlux48;
1380
1381   typedef shuffle_order_engine<minstd_rand0, 256> knuth_b;
1382
1383   typedef minstd_rand0 default_random_engine;
1384
1385   /**
1386    * A standard interface to a platform-specific non-deterministic
1387    * random number generator (if any are available).
1388    */
1389   class random_device
1390   {
1391   public:
1392     /** The type of the generated random value. */
1393     typedef unsigned int result_type;
1394
1395     // constructors, destructors and member functions
1396
1397 #ifdef _GLIBCXX_USE_RANDOM_TR1
1398
1399     explicit
1400     random_device(const std::string& __token = "/dev/urandom")
1401     {
1402       if ((__token != "/dev/urandom" && __token != "/dev/random")
1403           || !(_M_file = std::fopen(__token.c_str(), "rb")))
1404         std::__throw_runtime_error(__N("random_device::"
1405                                        "random_device(const std::string&)"));
1406     }
1407
1408     ~random_device()
1409     { std::fclose(_M_file); }
1410
1411 #else
1412
1413     explicit
1414     random_device(const std::string& __token = "mt19937")
1415     : _M_mt(_M_strtoul(__token)) { }
1416
1417   private:
1418     static unsigned long
1419     _M_strtoul(const std::string& __str)
1420     {
1421       unsigned long __ret = 5489UL;
1422       if (__str != "mt19937")
1423         {
1424           const char* __nptr = __str.c_str();
1425           char* __endptr;
1426           __ret = std::strtoul(__nptr, &__endptr, 0);
1427           if (*__nptr == '\0' || *__endptr != '\0')
1428             std::__throw_runtime_error(__N("random_device::_M_strtoul"
1429                                            "(const std::string&)"));
1430         }
1431       return __ret;
1432     }
1433
1434   public:
1435
1436 #endif
1437
1438     result_type
1439     min() const
1440     { return std::numeric_limits<result_type>::min(); }
1441
1442     result_type
1443     max() const
1444     { return std::numeric_limits<result_type>::max(); }
1445
1446     double
1447     entropy() const
1448     { return 0.0; }
1449
1450     result_type
1451     operator()()
1452     {
1453 #ifdef _GLIBCXX_USE_RANDOM_TR1
1454       result_type __ret;
1455       std::fread(reinterpret_cast<void*>(&__ret), sizeof(result_type),
1456                  1, _M_file);
1457       return __ret;
1458 #else
1459       return _M_mt();
1460 #endif
1461     }
1462
1463     // No copy functions.
1464     random_device(const random_device&) = delete;
1465     void operator=(const random_device&) = delete;
1466
1467   private:
1468
1469 #ifdef _GLIBCXX_USE_RANDOM_TR1
1470     FILE*        _M_file;
1471 #else
1472     mt19937      _M_mt;
1473 #endif
1474   };
1475
1476   /* @} */ // group random_generators
1477
1478   /**
1479    * @addtogroup random_distributions Random Number Distributions
1480    * @ingroup random
1481    * @{
1482    */
1483
1484   /**
1485    * @addtogroup random_distributions_uniform Uniform Distributions
1486    * @ingroup random_distributions
1487    * @{
1488    */
1489
1490   /**
1491    * @brief Uniform discrete distribution for random numbers.
1492    * A discrete random distribution on the range @f$[min, max]@f$ with equal
1493    * probability throughout the range.
1494    */
1495   template<typename _IntType = int>
1496     class uniform_int_distribution
1497     {
1498       static_assert(std::is_integral<_IntType>::value,
1499                     "template argument not an integral type");
1500
1501     public:
1502       /** The type of the range of the distribution. */
1503       typedef _IntType result_type;
1504       /** Parameter type. */
1505       struct param_type
1506       {
1507         typedef uniform_int_distribution<_IntType> distribution_type;
1508
1509         explicit
1510         param_type(_IntType __a = 0,
1511                    _IntType __b = std::numeric_limits<_IntType>::max())
1512         : _M_a(__a), _M_b(__b)
1513         {
1514           _GLIBCXX_DEBUG_ASSERT(_M_a <= _M_b);
1515         }
1516
1517         result_type
1518         a() const
1519         { return _M_a; }
1520
1521         result_type
1522         b() const
1523         { return _M_b; }
1524
1525       private:
1526         _IntType _M_a;
1527         _IntType _M_b;
1528       };
1529
1530     public:
1531       /**
1532        * @brief Constructs a uniform distribution object.
1533        */
1534       explicit
1535       uniform_int_distribution(_IntType __a = 0,
1536                            _IntType __b = std::numeric_limits<_IntType>::max())
1537       : _M_param(__a, __b)
1538       { }
1539
1540       explicit
1541       uniform_int_distribution(const param_type& __p)
1542       : _M_param(__p)
1543       { }
1544
1545       /**
1546        * @brief Resets the distribution state.
1547        *
1548        * Does nothing for the uniform integer distribution.
1549        */
1550       void
1551       reset() { }
1552
1553       result_type
1554       a() const
1555       { return _M_param.a(); }
1556
1557       result_type
1558       b() const
1559       { return _M_param.b(); }
1560
1561       /**
1562        * @brief Returns the inclusive lower bound of the distribution range.
1563        */
1564       result_type
1565       min() const
1566       { return this->a(); }
1567
1568       /**
1569        * @brief Returns the inclusive upper bound of the distribution range.
1570        */
1571       result_type
1572       max() const
1573       { return this->b(); }
1574
1575       /**
1576        * @brief Returns the parameter set of the distribution.
1577        */
1578       param_type
1579       param() const
1580       { return _M_param; }
1581
1582       /**
1583        * @brief Sets the parameter set of the distribution.
1584        * @param __param The new parameter set of the distribution.
1585        */
1586       void
1587       param(const param_type& __param)
1588       { _M_param = __param; }
1589
1590       /**
1591        * Gets a uniformly distributed random number in the range
1592        * @f$(min, max)@f$.
1593        */
1594       template<typename _UniformRandomNumberGenerator>
1595         result_type
1596         operator()(_UniformRandomNumberGenerator& __urng)
1597         { return this->operator()(__urng, this->param()); }
1598
1599       /**
1600        * Gets a uniform random number in the range @f$[0, n)@f$.
1601        *
1602        * This function is aimed at use with std::random_shuffle.
1603        */
1604       template<typename _UniformRandomNumberGenerator>
1605         result_type
1606         operator()(_UniformRandomNumberGenerator& __urng,
1607                    const param_type& __p);
1608
1609       param_type _M_param;
1610     };
1611
1612   /**
1613    * @brief Inserts a %uniform_int_distribution random number
1614    *        distribution @p __x into the output stream @p os.
1615    *
1616    * @param __os An output stream.
1617    * @param __x  A %uniform_int_distribution random number distribution.
1618    *
1619    * @returns The output stream with the state of @p __x inserted or in
1620    * an error state.
1621    */
1622   template<typename _IntType, typename _CharT, typename _Traits>
1623     std::basic_ostream<_CharT, _Traits>&
1624     operator<<(std::basic_ostream<_CharT, _Traits>&,
1625                const std::uniform_int_distribution<_IntType>&);
1626
1627   /**
1628    * @brief Extracts a %uniform_int_distribution random number distribution
1629    * @p __x from the input stream @p __is.
1630    *
1631    * @param __is An input stream.
1632    * @param __x  A %uniform_int_distribution random number generator engine.
1633    *
1634    * @returns The input stream with @p __x extracted or in an error state.
1635    */
1636   template<typename _IntType, typename _CharT, typename _Traits>
1637     std::basic_istream<_CharT, _Traits>&
1638     operator>>(std::basic_istream<_CharT, _Traits>&,
1639                std::uniform_int_distribution<_IntType>&);
1640
1641
1642   /**
1643    * @brief Uniform continuous distribution for random numbers.
1644    *
1645    * A continuous random distribution on the range [min, max) with equal
1646    * probability throughout the range.  The URNG should be real-valued and
1647    * deliver number in the range [0, 1).
1648    */
1649   template<typename _RealType = double>
1650     class uniform_real_distribution
1651     {
1652       static_assert(std::is_floating_point<_RealType>::value,
1653                     "template argument not a floating point type");
1654
1655     public:
1656       /** The type of the range of the distribution. */
1657       typedef _RealType result_type;
1658       /** Parameter type. */
1659       struct param_type
1660       {
1661         typedef uniform_real_distribution<_RealType> distribution_type;
1662
1663         explicit
1664         param_type(_RealType __a = _RealType(0),
1665                    _RealType __b = _RealType(1))
1666         : _M_a(__a), _M_b(__b)
1667         {
1668           _GLIBCXX_DEBUG_ASSERT(_M_a <= _M_b);
1669         }
1670
1671         result_type
1672         a() const
1673         { return _M_a; }
1674
1675         result_type
1676         b() const
1677         { return _M_b; }
1678
1679       private:
1680         _RealType _M_a;
1681         _RealType _M_b;
1682       };
1683
1684     public:
1685       /**
1686        * @brief Constructs a uniform_real_distribution object.
1687        *
1688        * @param __min [IN]  The lower bound of the distribution.
1689        * @param __max [IN]  The upper bound of the distribution.
1690        */
1691       explicit
1692       uniform_real_distribution(_RealType __a = _RealType(0),
1693                                 _RealType __b = _RealType(1))
1694       : _M_param(__a, __b)
1695       { }
1696
1697       explicit
1698       uniform_real_distribution(const param_type& __p)
1699       : _M_param(__p)
1700       { }
1701
1702       /**
1703        * @brief Resets the distribution state.
1704        *
1705        * Does nothing for the uniform real distribution.
1706        */
1707       void
1708       reset() { }
1709
1710       result_type
1711       a() const
1712       { return _M_param.a(); }
1713
1714       result_type
1715       b() const
1716       { return _M_param.b(); }
1717
1718       /**
1719        * @brief Returns the inclusive lower bound of the distribution range.
1720        */
1721       result_type
1722       min() const
1723       { return this->a(); }
1724
1725       /**
1726        * @brief Returns the inclusive upper bound of the distribution range.
1727        */
1728       result_type
1729       max() const
1730       { return this->b(); }
1731
1732       /**
1733        * @brief Returns the parameter set of the distribution.
1734        */
1735       param_type
1736       param() const
1737       { return _M_param; }
1738
1739       /**
1740        * @brief Sets the parameter set of the distribution.
1741        * @param __param The new parameter set of the distribution.
1742        */
1743       void
1744       param(const param_type& __param)
1745       { _M_param = __param; }
1746
1747       template<typename _UniformRandomNumberGenerator>
1748         result_type
1749         operator()(_UniformRandomNumberGenerator& __urng)
1750         { return this->operator()(__urng, this->param()); }
1751
1752       template<typename _UniformRandomNumberGenerator>
1753         result_type
1754         operator()(_UniformRandomNumberGenerator& __urng,
1755                    const param_type& __p)
1756         {
1757           __detail::_Adaptor<_UniformRandomNumberGenerator, result_type>
1758             __aurng(__urng);
1759           return (__aurng() * (__p.b() - __p.a())) + __p.a();
1760         }
1761
1762     private:
1763       param_type _M_param;
1764     };
1765
1766   /**
1767    * @brief Inserts a %uniform_real_distribution random number
1768    *        distribution @p __x into the output stream @p __os.
1769    *
1770    * @param __os An output stream.
1771    * @param __x  A %uniform_real_distribution random number distribution.
1772    *
1773    * @returns The output stream with the state of @p __x inserted or in
1774    *          an error state.
1775    */
1776   template<typename _RealType, typename _CharT, typename _Traits>
1777     std::basic_ostream<_CharT, _Traits>&
1778     operator<<(std::basic_ostream<_CharT, _Traits>&,
1779                const std::uniform_real_distribution<_RealType>&);
1780
1781   /**
1782    * @brief Extracts a %uniform_real_distribution random number distribution
1783    * @p __x from the input stream @p __is.
1784    *
1785    * @param __is An input stream.
1786    * @param __x  A %uniform_real_distribution random number generator engine.
1787    *
1788    * @returns The input stream with @p __x extracted or in an error state.
1789    */
1790   template<typename _RealType, typename _CharT, typename _Traits>
1791     std::basic_istream<_CharT, _Traits>&
1792     operator>>(std::basic_istream<_CharT, _Traits>&,
1793                std::uniform_real_distribution<_RealType>&);
1794
1795   /* @} */ // group random_distributions_uniform
1796
1797   /**
1798    * @addtogroup random_distributions_normal Normal Distributions
1799    * @ingroup random_distributions
1800    * @{
1801    */
1802
1803   /**
1804    * @brief A normal continuous distribution for random numbers.
1805    *
1806    * The formula for the normal probability density function is
1807    * @f[
1808    *     p(x|\mu,\sigma) = \frac{1}{\sigma \sqrt{2 \pi}}
1809    *            e^{- \frac{{x - \mu}^ {2}}{2 \sigma ^ {2}} }
1810    * @f]
1811    */
1812   template<typename _RealType = double>
1813     class normal_distribution
1814     {
1815       static_assert(std::is_floating_point<_RealType>::value,
1816                     "template argument not a floating point type");
1817
1818     public:
1819       /** The type of the range of the distribution. */
1820       typedef _RealType result_type;
1821       /** Parameter type. */
1822       struct param_type
1823       {
1824         typedef normal_distribution<_RealType> distribution_type;
1825
1826         explicit
1827         param_type(_RealType __mean = _RealType(0),
1828                    _RealType __stddev = _RealType(1))
1829         : _M_mean(__mean), _M_stddev(__stddev)
1830         {
1831           _GLIBCXX_DEBUG_ASSERT(_M_stddev > _RealType(0));
1832         }
1833
1834         _RealType
1835         mean() const
1836         { return _M_mean; }
1837
1838         _RealType
1839         stddev() const
1840         { return _M_stddev; }
1841
1842       private:
1843         _RealType _M_mean;
1844         _RealType _M_stddev;
1845       };
1846
1847     public:
1848       /**
1849        * Constructs a normal distribution with parameters @f$mean@f$ and
1850        * standard deviation.
1851        */
1852       explicit
1853       normal_distribution(result_type __mean = result_type(0),
1854                           result_type __stddev = result_type(1))
1855       : _M_param(__mean, __stddev), _M_saved_available(false)
1856       { }
1857
1858       explicit
1859       normal_distribution(const param_type& __p)
1860       : _M_param(__p), _M_saved_available(false)
1861       { }
1862
1863       /**
1864        * @brief Resets the distribution state.
1865        */
1866       void
1867       reset()
1868       { _M_saved_available = false; }
1869
1870       /**
1871        * @brief Returns the mean of the distribution.
1872        */
1873       _RealType
1874       mean() const
1875       { return _M_param.mean(); }
1876
1877       /**
1878        * @brief Returns the standard deviation of the distribution.
1879        */
1880       _RealType
1881       stddev() const
1882       { return _M_param.stddev(); }
1883
1884       /**
1885        * @brief Returns the parameter set of the distribution.
1886        */
1887       param_type
1888       param() const
1889       { return _M_param; }
1890
1891       /**
1892        * @brief Sets the parameter set of the distribution.
1893        * @param __param The new parameter set of the distribution.
1894        */
1895       void
1896       param(const param_type& __param)
1897       { _M_param = __param; }
1898
1899       /**
1900        * @brief Returns the greatest lower bound value of the distribution.
1901        */
1902       result_type
1903       min() const
1904       { return std::numeric_limits<result_type>::min(); }
1905
1906       /**
1907        * @brief Returns the least upper bound value of the distribution.
1908        */
1909       result_type
1910       max() const
1911       { return std::numeric_limits<result_type>::max(); }
1912
1913       template<typename _UniformRandomNumberGenerator>
1914         result_type
1915         operator()(_UniformRandomNumberGenerator& __urng)
1916         { return this->operator()(__urng, this->param()); }
1917
1918       template<typename _UniformRandomNumberGenerator>
1919         result_type
1920         operator()(_UniformRandomNumberGenerator& __urng,
1921                    const param_type& __p);
1922
1923       /**
1924        * @brief Inserts a %normal_distribution random number distribution
1925        * @p __x into the output stream @p __os.
1926        *
1927        * @param __os An output stream.
1928        * @param __x  A %normal_distribution random number distribution.
1929        *
1930        * @returns The output stream with the state of @p __x inserted or in
1931        * an error state.
1932        */
1933       template<typename _RealType1, typename _CharT, typename _Traits>
1934         friend std::basic_ostream<_CharT, _Traits>&
1935         operator<<(std::basic_ostream<_CharT, _Traits>&,
1936                    const std::normal_distribution<_RealType1>&);
1937
1938       /**
1939        * @brief Extracts a %normal_distribution random number distribution
1940        * @p __x from the input stream @p __is.
1941        *
1942        * @param __is An input stream.
1943        * @param __x  A %normal_distribution random number generator engine.
1944        *
1945        * @returns The input stream with @p __x extracted or in an error
1946        *          state.
1947        */
1948       template<typename _RealType1, typename _CharT, typename _Traits>
1949         friend std::basic_istream<_CharT, _Traits>&
1950         operator>>(std::basic_istream<_CharT, _Traits>&,
1951                    std::normal_distribution<_RealType1>&);
1952
1953     private:
1954       param_type  _M_param;
1955       result_type _M_saved;
1956       bool        _M_saved_available;
1957     };
1958
1959
1960   /**
1961    * @brief A lognormal_distribution random number distribution.
1962    *
1963    * The formula for the normal probability mass function is
1964    * @f[
1965    *     p(x|m,s) = \frac{1}{sx\sqrt{2\pi}}
1966    *                \exp{-\frac{(\ln{x} - m)^2}{2s^2}}
1967    * @f]
1968    */
1969   template<typename _RealType = double>
1970     class lognormal_distribution
1971     {
1972       static_assert(std::is_floating_point<_RealType>::value,
1973                     "template argument not a floating point type");
1974
1975     public:
1976       /** The type of the range of the distribution. */
1977       typedef _RealType result_type;
1978       /** Parameter type. */
1979       struct param_type
1980       {
1981         typedef lognormal_distribution<_RealType> distribution_type;
1982
1983         explicit
1984         param_type(_RealType __m = _RealType(0),
1985                    _RealType __s = _RealType(1))
1986         : _M_m(__m), _M_s(__s)
1987         { }
1988
1989         _RealType
1990         m() const
1991         { return _M_m; }
1992
1993         _RealType
1994         s() const
1995         { return _M_s; }
1996
1997       private:
1998         _RealType _M_m;
1999         _RealType _M_s;
2000       };
2001
2002       explicit
2003       lognormal_distribution(_RealType __m = _RealType(0),
2004                              _RealType __s = _RealType(1))
2005       : _M_param(__m, __s), _M_nd()
2006       { }
2007
2008       explicit
2009       lognormal_distribution(const param_type& __p)
2010       : _M_param(__p), _M_nd()
2011       { }
2012
2013       /**
2014        * Resets the distribution state.
2015        */
2016       void
2017       reset()
2018       { _M_nd.reset(); }
2019
2020       /**
2021        *
2022        */
2023       _RealType
2024       m() const
2025       { return _M_param.m(); }
2026
2027       _RealType
2028       s() const
2029       { return _M_param.s(); }
2030
2031       /**
2032        * @brief Returns the parameter set of the distribution.
2033        */
2034       param_type
2035       param() const
2036       { return _M_param; }
2037
2038       /**
2039        * @brief Sets the parameter set of the distribution.
2040        * @param __param The new parameter set of the distribution.
2041        */
2042       void
2043       param(const param_type& __param)
2044       { _M_param = __param; }
2045
2046       /**
2047        * @brief Returns the greatest lower bound value of the distribution.
2048        */
2049       result_type
2050       min() const
2051       { return result_type(0); }
2052
2053       /**
2054        * @brief Returns the least upper bound value of the distribution.
2055        */
2056       result_type
2057       max() const
2058       { return std::numeric_limits<result_type>::max(); }
2059
2060       template<typename _UniformRandomNumberGenerator>
2061         result_type
2062         operator()(_UniformRandomNumberGenerator& __urng)
2063         { return this->operator()(__urng, this->param()); }
2064
2065       template<typename _UniformRandomNumberGenerator>
2066         result_type
2067         operator()(_UniformRandomNumberGenerator& __urng,
2068                    const param_type& __p)
2069         { return std::exp(__p.s() * _M_nd(__urng) + __p.m()); }
2070
2071       /**
2072        * @brief Inserts a %lognormal_distribution random number distribution
2073        * @p __x into the output stream @p __os.
2074        *
2075        * @param __os An output stream.
2076        * @param __x  A %lognormal_distribution random number distribution.
2077        *
2078        * @returns The output stream with the state of @p __x inserted or in
2079        * an error state.
2080        */
2081       template<typename _RealType1, typename _CharT, typename _Traits>
2082         friend std::basic_ostream<_CharT, _Traits>&
2083         operator<<(std::basic_ostream<_CharT, _Traits>&,
2084                    const std::lognormal_distribution<_RealType1>&);
2085
2086       /**
2087        * @brief Extracts a %lognormal_distribution random number distribution
2088        * @p __x from the input stream @p __is.
2089        *
2090        * @param __is An input stream.
2091        * @param __x A %lognormal_distribution random number
2092        *            generator engine.
2093        *
2094        * @returns The input stream with @p __x extracted or in an error state.
2095        */
2096       template<typename _RealType1, typename _CharT, typename _Traits>
2097         friend std::basic_istream<_CharT, _Traits>&
2098         operator>>(std::basic_istream<_CharT, _Traits>&,
2099                    std::lognormal_distribution<_RealType1>&);
2100
2101     private:
2102       param_type _M_param;
2103
2104       std::normal_distribution<result_type> _M_nd;
2105     };
2106
2107
2108   /**
2109    * @brief A gamma continuous distribution for random numbers.
2110    *
2111    * The formula for the gamma probability density function is:
2112    * @f[
2113    *     p(x|\alpha,\beta) = \frac{1}{\beta\Gamma(\alpha)}
2114    *                         (x/\beta)^{\alpha - 1} e^{-x/\beta}
2115    * @f]
2116    */
2117   template<typename _RealType = double>
2118     class gamma_distribution
2119     {
2120       static_assert(std::is_floating_point<_RealType>::value,
2121                     "template argument not a floating point type");
2122
2123     public:
2124       /** The type of the range of the distribution. */
2125       typedef _RealType result_type;
2126       /** Parameter type. */
2127       struct param_type
2128       {
2129         typedef gamma_distribution<_RealType> distribution_type;
2130         friend class gamma_distribution<_RealType>;
2131
2132         explicit
2133         param_type(_RealType __alpha_val = _RealType(1),
2134                    _RealType __beta_val = _RealType(1))
2135         : _M_alpha(__alpha_val), _M_beta(__beta_val)
2136         {
2137           _GLIBCXX_DEBUG_ASSERT(_M_alpha > _RealType(0));
2138           _M_initialize();
2139         }
2140
2141         _RealType
2142         alpha() const
2143         { return _M_alpha; }
2144
2145         _RealType
2146         beta() const
2147         { return _M_beta; }
2148
2149       private:
2150         void
2151         _M_initialize();
2152
2153         _RealType _M_alpha;
2154         _RealType _M_beta;
2155
2156         _RealType _M_malpha, _M_a2;
2157       };
2158
2159     public:
2160       /**
2161        * @brief Constructs a gamma distribution with parameters
2162        * @f$\alpha@f$ and @f$\beta@f$.
2163        */
2164       explicit
2165       gamma_distribution(_RealType __alpha_val = _RealType(1),
2166                          _RealType __beta_val = _RealType(1))
2167       : _M_param(__alpha_val, __beta_val), _M_nd()
2168       { }
2169
2170       explicit
2171       gamma_distribution(const param_type& __p)
2172       : _M_param(__p), _M_nd()
2173       { }
2174
2175       /**
2176        * @brief Resets the distribution state.
2177        */
2178       void
2179       reset()
2180       { _M_nd.reset(); }
2181
2182       /**
2183        * @brief Returns the @f$\alpha@f$ of the distribution.
2184        */
2185       _RealType
2186       alpha() const
2187       { return _M_param.alpha(); }
2188
2189       /**
2190        * @brief Returns the @f$\beta@f$ of the distribution.
2191        */
2192       _RealType
2193       beta() const
2194       { return _M_param.beta(); }
2195
2196       /**
2197        * @brief Returns the parameter set of the distribution.
2198        */
2199       param_type
2200       param() const
2201       { return _M_param; }
2202
2203       /**
2204        * @brief Sets the parameter set of the distribution.
2205        * @param __param The new parameter set of the distribution.
2206        */
2207       void
2208       param(const param_type& __param)
2209       { _M_param = __param; }
2210
2211       /**
2212        * @brief Returns the greatest lower bound value of the distribution.
2213        */
2214       result_type
2215       min() const
2216       { return result_type(0); }
2217
2218       /**
2219        * @brief Returns the least upper bound value of the distribution.
2220        */
2221       result_type
2222       max() const
2223       { return std::numeric_limits<result_type>::max(); }
2224
2225       template<typename _UniformRandomNumberGenerator>
2226         result_type
2227         operator()(_UniformRandomNumberGenerator& __urng)
2228         { return this->operator()(__urng, this->param()); }
2229
2230       template<typename _UniformRandomNumberGenerator>
2231         result_type
2232         operator()(_UniformRandomNumberGenerator& __urng,
2233                    const param_type& __p);
2234
2235       /**
2236        * @brief Inserts a %gamma_distribution random number distribution
2237        * @p __x into the output stream @p __os.
2238        *
2239        * @param __os An output stream.
2240        * @param __x  A %gamma_distribution random number distribution.
2241        *
2242        * @returns The output stream with the state of @p __x inserted or in
2243        * an error state.
2244        */
2245       template<typename _RealType1, typename _CharT, typename _Traits>
2246         friend std::basic_ostream<_CharT, _Traits>&
2247         operator<<(std::basic_ostream<_CharT, _Traits>&,
2248                    const std::gamma_distribution<_RealType1>&);
2249
2250       /**
2251        * @brief Extracts a %gamma_distribution random number distribution
2252        * @p __x from the input stream @p __is.
2253        *
2254        * @param __is An input stream.
2255        * @param __x  A %gamma_distribution random number generator engine.
2256        *
2257        * @returns The input stream with @p __x extracted or in an error state.
2258        */
2259       template<typename _RealType1, typename _CharT, typename _Traits>
2260         friend std::basic_istream<_CharT, _Traits>&
2261         operator>>(std::basic_istream<_CharT, _Traits>&,
2262                    std::gamma_distribution<_RealType1>&);
2263
2264     private:
2265       param_type _M_param;
2266
2267       std::normal_distribution<result_type> _M_nd;
2268     };
2269
2270
2271   /**
2272    * @brief A chi_squared_distribution random number distribution.
2273    *
2274    * The formula for the normal probability mass function is
2275    * @f$p(x|n) = \frac{x^{(n/2) - 1}e^{-x/2}}{\Gamma(n/2) 2^{n/2}}@f$
2276    */
2277   template<typename _RealType = double>
2278     class chi_squared_distribution
2279     {
2280       static_assert(std::is_floating_point<_RealType>::value,
2281                     "template argument not a floating point type");
2282
2283     public:
2284       /** The type of the range of the distribution. */
2285       typedef _RealType result_type;
2286       /** Parameter type. */
2287       struct param_type
2288       {
2289         typedef chi_squared_distribution<_RealType> distribution_type;
2290
2291         explicit
2292         param_type(_RealType __n = _RealType(1))
2293         : _M_n(__n)
2294         { }
2295
2296         _RealType
2297         n() const
2298         { return _M_n; }
2299
2300       private:
2301         _RealType _M_n;
2302       };
2303
2304       explicit
2305       chi_squared_distribution(_RealType __n = _RealType(1))
2306       : _M_param(__n), _M_gd(__n / 2)
2307       { }
2308
2309       explicit
2310       chi_squared_distribution(const param_type& __p)
2311       : _M_param(__p), _M_gd(__p.n() / 2)
2312       { }
2313
2314       /**
2315        * @brief Resets the distribution state.
2316        */
2317       void
2318       reset()
2319       { _M_gd.reset(); }
2320
2321       /**
2322        *
2323        */
2324       _RealType
2325       n() const
2326       { return _M_param.n(); }
2327
2328       /**
2329        * @brief Returns the parameter set of the distribution.
2330        */
2331       param_type
2332       param() const
2333       { return _M_param; }
2334
2335       /**
2336        * @brief Sets the parameter set of the distribution.
2337        * @param __param The new parameter set of the distribution.
2338        */
2339       void
2340       param(const param_type& __param)
2341       { _M_param = __param; }
2342
2343       /**
2344        * @brief Returns the greatest lower bound value of the distribution.
2345        */
2346       result_type
2347       min() const
2348       { return result_type(0); }
2349
2350       /**
2351        * @brief Returns the least upper bound value of the distribution.
2352        */
2353       result_type
2354       max() const
2355       { return std::numeric_limits<result_type>::max(); }
2356
2357       template<typename _UniformRandomNumberGenerator>
2358         result_type
2359         operator()(_UniformRandomNumberGenerator& __urng)
2360         { return 2 * _M_gd(__urng); }
2361
2362       template<typename _UniformRandomNumberGenerator>
2363         result_type
2364         operator()(_UniformRandomNumberGenerator& __urng,
2365                    const param_type& __p)
2366         {
2367           typedef typename std::gamma_distribution<result_type>::param_type
2368             param_type;
2369           return 2 * _M_gd(__urng, param_type(__p.n() / 2));
2370         }
2371
2372       /**
2373        * @brief Inserts a %chi_squared_distribution random number distribution
2374        * @p __x into the output stream @p __os.
2375        *
2376        * @param __os An output stream.
2377        * @param __x  A %chi_squared_distribution random number distribution.
2378        *
2379        * @returns The output stream with the state of @p __x inserted or in
2380        * an error state.
2381        */
2382       template<typename _RealType1, typename _CharT, typename _Traits>
2383         friend std::basic_ostream<_CharT, _Traits>&
2384         operator<<(std::basic_ostream<_CharT, _Traits>&,
2385                    const std::chi_squared_distribution<_RealType1>&);
2386
2387       /**
2388        * @brief Extracts a %chi_squared_distribution random number distribution
2389        * @p __x from the input stream @p __is.
2390        *
2391        * @param __is An input stream.
2392        * @param __x A %chi_squared_distribution random number
2393        *            generator engine.
2394        *
2395        * @returns The input stream with @p __x extracted or in an error state.
2396        */
2397       template<typename _RealType1, typename _CharT, typename _Traits>
2398         friend std::basic_istream<_CharT, _Traits>&
2399         operator>>(std::basic_istream<_CharT, _Traits>&,
2400                    std::chi_squared_distribution<_RealType1>&);
2401
2402     private:
2403       param_type _M_param;
2404
2405       std::gamma_distribution<result_type> _M_gd;
2406     };
2407
2408
2409   /**
2410    * @brief A cauchy_distribution random number distribution.
2411    *
2412    * The formula for the normal probability mass function is
2413    * @f$p(x|a,b) = (\pi b (1 + (\frac{x-a}{b})^2))^{-1}@f$
2414    */
2415   template<typename _RealType = double>
2416     class cauchy_distribution
2417     {
2418       static_assert(std::is_floating_point<_RealType>::value,
2419                     "template argument not a floating point type");
2420
2421     public:
2422       /** The type of the range of the distribution. */
2423       typedef _RealType result_type;
2424       /** Parameter type. */
2425       struct param_type
2426       {
2427         typedef cauchy_distribution<_RealType> distribution_type;
2428
2429         explicit
2430         param_type(_RealType __a = _RealType(0),
2431                    _RealType __b = _RealType(1))
2432         : _M_a(__a), _M_b(__b)
2433         { }
2434
2435         _RealType
2436         a() const
2437         { return _M_a; }
2438
2439         _RealType
2440         b() const
2441         { return _M_b; }
2442
2443       private:
2444         _RealType _M_a;
2445         _RealType _M_b;
2446       };
2447
2448       explicit
2449       cauchy_distribution(_RealType __a = _RealType(0),
2450                           _RealType __b = _RealType(1))
2451       : _M_param(__a, __b)
2452       { }
2453
2454       explicit
2455       cauchy_distribution(const param_type& __p)
2456       : _M_param(__p)
2457       { }
2458
2459       /**
2460        * @brief Resets the distribution state.
2461        */
2462       void
2463       reset()
2464       { }
2465
2466       /**
2467        *
2468        */
2469       _RealType
2470       a() const
2471       { return _M_param.a(); }
2472
2473       _RealType
2474       b() const
2475       { return _M_param.b(); }
2476
2477       /**
2478        * @brief Returns the parameter set of the distribution.
2479        */
2480       param_type
2481       param() const
2482       { return _M_param; }
2483
2484       /**
2485        * @brief Sets the parameter set of the distribution.
2486        * @param __param The new parameter set of the distribution.
2487        */
2488       void
2489       param(const param_type& __param)
2490       { _M_param = __param; }
2491
2492       /**
2493        * @brief Returns the greatest lower bound value of the distribution.
2494        */
2495       result_type
2496       min() const
2497       { return std::numeric_limits<result_type>::min(); }
2498
2499       /**
2500        * @brief Returns the least upper bound value of the distribution.
2501        */
2502       result_type
2503       max() const
2504       { return std::numeric_limits<result_type>::max(); }
2505
2506       template<typename _UniformRandomNumberGenerator>
2507         result_type
2508         operator()(_UniformRandomNumberGenerator& __urng)
2509         { return this->operator()(__urng, this->param()); }
2510
2511       template<typename _UniformRandomNumberGenerator>
2512         result_type
2513         operator()(_UniformRandomNumberGenerator& __urng,
2514                    const param_type& __p);
2515
2516     private:
2517       param_type _M_param;
2518     };
2519
2520   /**
2521    * @brief Inserts a %cauchy_distribution random number distribution
2522    * @p __x into the output stream @p __os.
2523    *
2524    * @param __os An output stream.
2525    * @param __x  A %cauchy_distribution random number distribution.
2526    *
2527    * @returns The output stream with the state of @p __x inserted or in
2528    * an error state.
2529    */
2530   template<typename _RealType, typename _CharT, typename _Traits>
2531     std::basic_ostream<_CharT, _Traits>&
2532     operator<<(std::basic_ostream<_CharT, _Traits>&,
2533                const std::cauchy_distribution<_RealType>&);
2534
2535   /**
2536    * @brief Extracts a %cauchy_distribution random number distribution
2537    * @p __x from the input stream @p __is.
2538    *
2539    * @param __is An input stream.
2540    * @param __x A %cauchy_distribution random number
2541    *            generator engine.
2542    *
2543    * @returns The input stream with @p __x extracted or in an error state.
2544    */
2545   template<typename _RealType, typename _CharT, typename _Traits>
2546     std::basic_istream<_CharT, _Traits>&
2547     operator>>(std::basic_istream<_CharT, _Traits>&,
2548                std::cauchy_distribution<_RealType>&);
2549
2550
2551   /**
2552    * @brief A fisher_f_distribution random number distribution.
2553    *
2554    * The formula for the normal probability mass function is
2555    * @f[
2556    *     p(x|m,n) = \frac{\Gamma((m+n)/2)}{\Gamma(m/2)\Gamma(n/2)}
2557    *                (\frac{m}{n})^{m/2} x^{(m/2)-1}
2558    *                (1 + \frac{mx}{n})^{-(m+n)/2}
2559    * @f]
2560    */
2561   template<typename _RealType = double>
2562     class fisher_f_distribution
2563     {
2564       static_assert(std::is_floating_point<_RealType>::value,
2565                     "template argument not a floating point type");
2566
2567     public:
2568       /** The type of the range of the distribution. */
2569       typedef _RealType result_type;
2570       /** Parameter type. */
2571       struct param_type
2572       {
2573         typedef fisher_f_distribution<_RealType> distribution_type;
2574
2575         explicit
2576         param_type(_RealType __m = _RealType(1),
2577                    _RealType __n = _RealType(1))
2578         : _M_m(__m), _M_n(__n)
2579         { }
2580
2581         _RealType
2582         m() const
2583         { return _M_m; }
2584
2585         _RealType
2586         n() const
2587         { return _M_n; }
2588
2589       private:
2590         _RealType _M_m;
2591         _RealType _M_n;
2592       };
2593
2594       explicit
2595       fisher_f_distribution(_RealType __m = _RealType(1),
2596                             _RealType __n = _RealType(1))
2597       : _M_param(__m, __n), _M_gd_x(__m / 2), _M_gd_y(__n / 2)
2598       { }
2599
2600       explicit
2601       fisher_f_distribution(const param_type& __p)
2602       : _M_param(__p), _M_gd_x(__p.m() / 2), _M_gd_y(__p.n() / 2)
2603       { }
2604
2605       /**
2606        * @brief Resets the distribution state.
2607        */
2608       void
2609       reset()
2610       {
2611         _M_gd_x.reset();
2612         _M_gd_y.reset();
2613       }
2614
2615       /**
2616        *
2617        */
2618       _RealType
2619       m() const
2620       { return _M_param.m(); }
2621
2622       _RealType
2623       n() const
2624       { return _M_param.n(); }
2625
2626       /**
2627        * @brief Returns the parameter set of the distribution.
2628        */
2629       param_type
2630       param() const
2631       { return _M_param; }
2632
2633       /**
2634        * @brief Sets the parameter set of the distribution.
2635        * @param __param The new parameter set of the distribution.
2636        */
2637       void
2638       param(const param_type& __param)
2639       { _M_param = __param; }
2640
2641       /**
2642        * @brief Returns the greatest lower bound value of the distribution.
2643        */
2644       result_type
2645       min() const
2646       { return result_type(0); }
2647
2648       /**
2649        * @brief Returns the least upper bound value of the distribution.
2650        */
2651       result_type
2652       max() const
2653       { return std::numeric_limits<result_type>::max(); }
2654
2655       template<typename _UniformRandomNumberGenerator>
2656         result_type
2657         operator()(_UniformRandomNumberGenerator& __urng)
2658         { return (_M_gd_x(__urng) * n()) / (_M_gd_y(__urng) * m()); }
2659
2660       template<typename _UniformRandomNumberGenerator>
2661         result_type
2662         operator()(_UniformRandomNumberGenerator& __urng,
2663                    const param_type& __p)
2664         {
2665           typedef typename std::gamma_distribution<result_type>::param_type
2666             param_type;
2667           return ((_M_gd_x(__urng, param_type(__p.m() / 2)) * n())
2668                   / (_M_gd_y(__urng, param_type(__p.n() / 2)) * m()));
2669         }
2670
2671       /**
2672        * @brief Inserts a %fisher_f_distribution random number distribution
2673        * @p __x into the output stream @p __os.
2674        *
2675        * @param __os An output stream.
2676        * @param __x  A %fisher_f_distribution random number distribution.
2677        *
2678        * @returns The output stream with the state of @p __x inserted or in
2679        * an error state.
2680        */
2681       template<typename _RealType1, typename _CharT, typename _Traits>
2682         friend std::basic_ostream<_CharT, _Traits>&
2683         operator<<(std::basic_ostream<_CharT, _Traits>&,
2684                    const std::fisher_f_distribution<_RealType1>&);
2685
2686       /**
2687        * @brief Extracts a %fisher_f_distribution random number distribution
2688        * @p __x from the input stream @p __is.
2689        *
2690        * @param __is An input stream.
2691        * @param __x A %fisher_f_distribution random number
2692        *            generator engine.
2693        *
2694        * @returns The input stream with @p __x extracted or in an error state.
2695        */
2696       template<typename _RealType1, typename _CharT, typename _Traits>
2697         friend std::basic_istream<_CharT, _Traits>&
2698         operator>>(std::basic_istream<_CharT, _Traits>&,
2699                    std::fisher_f_distribution<_RealType1>&);
2700
2701     private:
2702       param_type _M_param;
2703
2704       std::gamma_distribution<result_type> _M_gd_x, _M_gd_y;
2705     };
2706
2707
2708   /**
2709    * @brief A student_t_distribution random number distribution.
2710    *
2711    * The formula for the normal probability mass function is:
2712    * @f[
2713    *     p(x|n) = \frac{1}{\sqrt(n\pi)} \frac{\Gamma((n+1)/2)}{\Gamma(n/2)}
2714    *              (1 + \frac{x^2}{n}) ^{-(n+1)/2}
2715    * @f]
2716    */
2717   template<typename _RealType = double>
2718     class student_t_distribution
2719     {
2720       static_assert(std::is_floating_point<_RealType>::value,
2721                     "template argument not a floating point type");
2722
2723     public:
2724       /** The type of the range of the distribution. */
2725       typedef _RealType result_type;
2726       /** Parameter type. */
2727       struct param_type
2728       {
2729         typedef student_t_distribution<_RealType> distribution_type;
2730
2731         explicit
2732         param_type(_RealType __n = _RealType(1))
2733         : _M_n(__n)
2734         { }
2735
2736         _RealType
2737         n() const
2738         { return _M_n; }
2739
2740       private:
2741         _RealType _M_n;
2742       };
2743
2744       explicit
2745       student_t_distribution(_RealType __n = _RealType(1))
2746       : _M_param(__n), _M_nd(), _M_gd(__n / 2, 2)
2747       { }
2748
2749       explicit
2750       student_t_distribution(const param_type& __p)
2751       : _M_param(__p), _M_nd(), _M_gd(__p.n() / 2, 2)
2752       { }
2753
2754       /**
2755        * @brief Resets the distribution state.
2756        */
2757       void
2758       reset()
2759       {
2760         _M_nd.reset();
2761         _M_gd.reset();
2762       }
2763
2764       /**
2765        *
2766        */
2767       _RealType
2768       n() const
2769       { return _M_param.n(); }
2770
2771       /**
2772        * @brief Returns the parameter set of the distribution.
2773        */
2774       param_type
2775       param() const
2776       { return _M_param; }
2777
2778       /**
2779        * @brief Sets the parameter set of the distribution.
2780        * @param __param The new parameter set of the distribution.
2781        */
2782       void
2783       param(const param_type& __param)
2784       { _M_param = __param; }
2785
2786       /**
2787        * @brief Returns the greatest lower bound value of the distribution.
2788        */
2789       result_type
2790       min() const
2791       { return std::numeric_limits<result_type>::min(); }
2792
2793       /**
2794        * @brief Returns the least upper bound value of the distribution.
2795        */
2796       result_type
2797       max() const
2798       { return std::numeric_limits<result_type>::max(); }
2799
2800       template<typename _UniformRandomNumberGenerator>
2801         result_type
2802         operator()(_UniformRandomNumberGenerator& __urng)
2803         { return _M_nd(__urng) * std::sqrt(n() / _M_gd(__urng)); }
2804
2805       template<typename _UniformRandomNumberGenerator>
2806         result_type
2807         operator()(_UniformRandomNumberGenerator& __urng,
2808                    const param_type& __p)
2809         {
2810           typedef typename std::gamma_distribution<result_type>::param_type
2811             param_type;
2812
2813           const result_type __g = _M_gd(__urng, param_type(__p.n() / 2, 2));
2814           return _M_nd(__urng) * std::sqrt(__p.n() / __g);
2815         }
2816
2817       /**
2818        * @brief Inserts a %student_t_distribution random number distribution
2819        * @p __x into the output stream @p __os.
2820        *
2821        * @param __os An output stream.
2822        * @param __x  A %student_t_distribution random number distribution.
2823        *
2824        * @returns The output stream with the state of @p __x inserted or in
2825        * an error state.
2826        */
2827       template<typename _RealType1, typename _CharT, typename _Traits>
2828         friend std::basic_ostream<_CharT, _Traits>&
2829         operator<<(std::basic_ostream<_CharT, _Traits>&,
2830                    const std::student_t_distribution<_RealType1>&);
2831
2832       /**
2833        * @brief Extracts a %student_t_distribution random number distribution
2834        * @p __x from the input stream @p __is.
2835        *
2836        * @param __is An input stream.
2837        * @param __x A %student_t_distribution random number
2838        *            generator engine.
2839        *
2840        * @returns The input stream with @p __x extracted or in an error state.
2841        */
2842       template<typename _RealType1, typename _CharT, typename _Traits>
2843         friend std::basic_istream<_CharT, _Traits>&
2844         operator>>(std::basic_istream<_CharT, _Traits>&,
2845                    std::student_t_distribution<_RealType1>&);
2846
2847     private:
2848       param_type _M_param;
2849
2850       std::normal_distribution<result_type> _M_nd;
2851       std::gamma_distribution<result_type> _M_gd;
2852     };
2853
2854   /* @} */ // group random_distributions_normal
2855
2856   /**
2857    * @addtogroup random_distributions_bernoulli Bernoulli Distributions
2858    * @ingroup random_distributions
2859    * @{
2860    */
2861
2862   /**
2863    * @brief A Bernoulli random number distribution.
2864    *
2865    * Generates a sequence of true and false values with likelihood @f$p@f$
2866    * that true will come up and @f$(1 - p)@f$ that false will appear.
2867    */
2868   class bernoulli_distribution
2869   {
2870   public:
2871     /** The type of the range of the distribution. */
2872     typedef bool result_type;
2873     /** Parameter type. */
2874     struct param_type
2875     {
2876       typedef bernoulli_distribution distribution_type;
2877
2878       explicit
2879       param_type(double __p = 0.5)
2880       : _M_p(__p)
2881       {
2882         _GLIBCXX_DEBUG_ASSERT((_M_p >= 0.0) && (_M_p <= 1.0));
2883       }
2884
2885       double
2886       p() const
2887       { return _M_p; }
2888
2889     private:
2890       double _M_p;
2891     };
2892
2893   public:
2894     /**
2895      * @brief Constructs a Bernoulli distribution with likelihood @p p.
2896      *
2897      * @param __p  [IN]  The likelihood of a true result being returned.
2898      *                   Must be in the interval @f$[0, 1]@f$.
2899      */
2900     explicit
2901     bernoulli_distribution(double __p = 0.5)
2902     : _M_param(__p)
2903     { }
2904
2905     explicit
2906     bernoulli_distribution(const param_type& __p)
2907     : _M_param(__p)
2908     { }
2909
2910     /**
2911      * @brief Resets the distribution state.
2912      *
2913      * Does nothing for a Bernoulli distribution.
2914      */
2915     void
2916     reset() { }
2917
2918     /**
2919      * @brief Returns the @p p parameter of the distribution.
2920      */
2921     double
2922     p() const
2923     { return _M_param.p(); }
2924
2925     /**
2926      * @brief Returns the parameter set of the distribution.
2927      */
2928     param_type
2929     param() const
2930     { return _M_param; }
2931
2932     /**
2933      * @brief Sets the parameter set of the distribution.
2934      * @param __param The new parameter set of the distribution.
2935      */
2936     void
2937     param(const param_type& __param)
2938     { _M_param = __param; }
2939
2940     /**
2941      * @brief Returns the greatest lower bound value of the distribution.
2942      */
2943     result_type
2944     min() const
2945     { return std::numeric_limits<result_type>::min(); }
2946
2947     /**
2948      * @brief Returns the least upper bound value of the distribution.
2949      */
2950     result_type
2951     max() const
2952     { return std::numeric_limits<result_type>::max(); }
2953
2954     /**
2955      * @brief Returns the next value in the Bernoullian sequence.
2956      */
2957     template<typename _UniformRandomNumberGenerator>
2958       result_type
2959       operator()(_UniformRandomNumberGenerator& __urng)
2960       { return this->operator()(__urng, this->param()); }
2961
2962     template<typename _UniformRandomNumberGenerator>
2963       result_type
2964       operator()(_UniformRandomNumberGenerator& __urng,
2965                  const param_type& __p)
2966       {
2967         __detail::_Adaptor<_UniformRandomNumberGenerator, double>
2968           __aurng(__urng);
2969         if ((__aurng() - __aurng.min())
2970              < __p.p() * (__aurng.max() - __aurng.min()))
2971           return true;
2972         return false;
2973       }
2974
2975   private:
2976     param_type _M_param;
2977   };
2978
2979   /**
2980    * @brief Inserts a %bernoulli_distribution random number distribution
2981    * @p __x into the output stream @p __os.
2982    *
2983    * @param __os An output stream.
2984    * @param __x  A %bernoulli_distribution random number distribution.
2985    *
2986    * @returns The output stream with the state of @p __x inserted or in
2987    * an error state.
2988    */
2989   template<typename _CharT, typename _Traits>
2990     std::basic_ostream<_CharT, _Traits>&
2991     operator<<(std::basic_ostream<_CharT, _Traits>&,
2992                const std::bernoulli_distribution&);
2993
2994   /**
2995    * @brief Extracts a %bernoulli_distribution random number distribution
2996    * @p __x from the input stream @p __is.
2997    *
2998    * @param __is An input stream.
2999    * @param __x  A %bernoulli_distribution random number generator engine.
3000    *
3001    * @returns The input stream with @p __x extracted or in an error state.
3002    */
3003   template<typename _CharT, typename _Traits>
3004     std::basic_istream<_CharT, _Traits>&
3005     operator>>(std::basic_istream<_CharT, _Traits>& __is,
3006                std::bernoulli_distribution& __x)
3007     {
3008       double __p;
3009       __is >> __p;
3010       __x.param(bernoulli_distribution::param_type(__p));
3011       return __is;
3012     }
3013
3014
3015   /**
3016    * @brief A discrete binomial random number distribution.
3017    *
3018    * The formula for the binomial probability density function is
3019    * @f$p(i|t,p) = \binom{n}{i} p^i (1 - p)^{t - i}@f$ where @f$t@f$
3020    * and @f$p@f$ are the parameters of the distribution.
3021    */
3022   template<typename _IntType = int>
3023     class binomial_distribution
3024     {
3025       static_assert(std::is_integral<_IntType>::value,
3026                     "template argument not an integral type");
3027
3028     public:
3029       /** The type of the range of the distribution. */
3030       typedef _IntType result_type;
3031       /** Parameter type. */
3032       struct param_type
3033       {
3034         typedef binomial_distribution<_IntType> distribution_type;
3035         friend class binomial_distribution<_IntType>;
3036
3037         explicit
3038         param_type(_IntType __t = _IntType(1), double __p = 0.5)
3039         : _M_t(__t), _M_p(__p)
3040         {
3041           _GLIBCXX_DEBUG_ASSERT((_M_t >= _IntType(0))
3042                                 && (_M_p >= 0.0)
3043                                 && (_M_p <= 1.0));
3044           _M_initialize();
3045         }
3046
3047         _IntType
3048         t() const
3049         { return _M_t; }
3050
3051         double
3052         p() const
3053         { return _M_p; }
3054
3055       private:
3056         void
3057         _M_initialize();
3058
3059         _IntType _M_t;
3060         double _M_p;
3061
3062         double _M_q;
3063 #if _GLIBCXX_USE_C99_MATH_TR1
3064         double _M_d1, _M_d2, _M_s1, _M_s2, _M_c,
3065                _M_a1, _M_a123, _M_s, _M_lf, _M_lp1p;
3066 #endif
3067         bool   _M_easy;
3068       };
3069
3070       // constructors and member function
3071       explicit
3072       binomial_distribution(_IntType __t = _IntType(1),
3073                             double __p = 0.5)
3074       : _M_param(__t, __p), _M_nd()
3075       { }
3076
3077       explicit
3078       binomial_distribution(const param_type& __p)
3079       : _M_param(__p), _M_nd()
3080       { }
3081
3082       /**
3083        * @brief Resets the distribution state.
3084        */
3085       void
3086       reset()
3087       { _M_nd.reset(); }
3088
3089       /**
3090        * @brief Returns the distribution @p t parameter.
3091        */
3092       _IntType
3093       t() const
3094       { return _M_param.t(); }
3095
3096       /**
3097        * @brief Returns the distribution @p p parameter.
3098        */
3099       double
3100       p() const
3101       { return _M_param.p(); }
3102
3103       /**
3104        * @brief Returns the parameter set of the distribution.
3105        */
3106       param_type
3107       param() const
3108       { return _M_param; }
3109
3110       /**
3111        * @brief Sets the parameter set of the distribution.
3112        * @param __param The new parameter set of the distribution.
3113        */
3114       void
3115       param(const param_type& __param)
3116       { _M_param = __param; }
3117
3118       /**
3119        * @brief Returns the greatest lower bound value of the distribution.
3120        */
3121       result_type
3122       min() const
3123       { return 0; }
3124
3125       /**
3126        * @brief Returns the least upper bound value of the distribution.
3127        */
3128       result_type
3129       max() const
3130       { return _M_param.t(); }
3131
3132       template<typename _UniformRandomNumberGenerator>
3133         result_type
3134         operator()(_UniformRandomNumberGenerator& __urng)
3135         { return this->operator()(__urng, this->param()); }
3136
3137       template<typename _UniformRandomNumberGenerator>
3138         result_type
3139         operator()(_UniformRandomNumberGenerator& __urng,
3140                    const param_type& __p);
3141
3142       /**
3143        * @brief Inserts a %binomial_distribution random number distribution
3144        * @p __x into the output stream @p __os.
3145        *
3146        * @param __os An output stream.
3147        * @param __x  A %binomial_distribution random number distribution.
3148        *
3149        * @returns The output stream with the state of @p __x inserted or in
3150        * an error state.
3151        */
3152       template<typename _IntType1,
3153                typename _CharT, typename _Traits>
3154         friend std::basic_ostream<_CharT, _Traits>&
3155         operator<<(std::basic_ostream<_CharT, _Traits>&,
3156                    const std::binomial_distribution<_IntType1>&);
3157
3158       /**
3159        * @brief Extracts a %binomial_distribution random number distribution
3160        * @p __x from the input stream @p __is.
3161        *
3162        * @param __is An input stream.
3163        * @param __x  A %binomial_distribution random number generator engine.
3164        *
3165        * @returns The input stream with @p __x extracted or in an error
3166        *          state.
3167        */
3168       template<typename _IntType1,
3169                typename _CharT, typename _Traits>
3170         friend std::basic_istream<_CharT, _Traits>&
3171         operator>>(std::basic_istream<_CharT, _Traits>&,
3172                    std::binomial_distribution<_IntType1>&);
3173
3174     private:
3175       template<typename _UniformRandomNumberGenerator>
3176         result_type
3177         _M_waiting(_UniformRandomNumberGenerator& __urng, _IntType __t);
3178
3179       param_type _M_param;
3180
3181       // NB: Unused when _GLIBCXX_USE_C99_MATH_TR1 is undefined.
3182       std::normal_distribution<double> _M_nd;
3183     };
3184
3185
3186   /**
3187    * @brief A discrete geometric random number distribution.
3188    *
3189    * The formula for the geometric probability density function is
3190    * @f$p(i|p) = (1 - p)p^{i-1}@f$ where @f$p@f$ is the parameter of the
3191    * distribution.
3192    */
3193   template<typename _IntType = int>
3194     class geometric_distribution
3195     {
3196       static_assert(std::is_integral<_IntType>::value,
3197                     "template argument not an integral type");
3198
3199     public:
3200       /** The type of the range of the distribution. */
3201       typedef _IntType  result_type;
3202       /** Parameter type. */
3203       struct param_type
3204       {
3205         typedef geometric_distribution<_IntType> distribution_type;
3206         friend class geometric_distribution<_IntType>;
3207
3208         explicit
3209         param_type(double __p = 0.5)
3210         : _M_p(__p)
3211         {
3212           _GLIBCXX_DEBUG_ASSERT((_M_p >= 0.0)
3213                              && (_M_p <= 1.0));
3214           _M_initialize();
3215         }
3216
3217         double
3218         p() const
3219         { return _M_p; }
3220
3221       private:
3222         void
3223         _M_initialize()
3224         { _M_log_p = std::log(_M_p); }
3225
3226         double _M_p;
3227
3228         double _M_log_p;
3229       };
3230
3231       // constructors and member function
3232       explicit
3233       geometric_distribution(double __p = 0.5)
3234       : _M_param(__p)
3235       { }
3236
3237       explicit
3238       geometric_distribution(const param_type& __p)
3239       : _M_param(__p)
3240       { }
3241
3242       /**
3243        * @brief Resets the distribution state.
3244        *
3245        * Does nothing for the geometric distribution.
3246        */
3247       void
3248       reset() { }
3249
3250       /**
3251        * @brief Returns the distribution parameter @p p.
3252        */
3253       double
3254       p() const
3255       { return _M_param.p(); }
3256
3257       /**
3258        * @brief Returns the parameter set of the distribution.
3259        */
3260       param_type
3261       param() const
3262       { return _M_param; }
3263
3264       /**
3265        * @brief Sets the parameter set of the distribution.
3266        * @param __param The new parameter set of the distribution.
3267        */
3268       void
3269       param(const param_type& __param)
3270       { _M_param = __param; }
3271
3272       /**
3273        * @brief Returns the greatest lower bound value of the distribution.
3274        */
3275       result_type
3276       min() const
3277       { return 0; }
3278
3279       /**
3280        * @brief Returns the least upper bound value of the distribution.
3281        */
3282       result_type
3283       max() const
3284       { return std::numeric_limits<result_type>::max(); }
3285
3286       template<typename _UniformRandomNumberGenerator>
3287         result_type
3288         operator()(_UniformRandomNumberGenerator& __urng)
3289         { return this->operator()(__urng, this->param()); }
3290
3291       template<typename _UniformRandomNumberGenerator>
3292         result_type
3293         operator()(_UniformRandomNumberGenerator& __urng,
3294                    const param_type& __p);
3295
3296     private:
3297       param_type _M_param;
3298     };
3299
3300   /**
3301    * @brief Inserts a %geometric_distribution random number distribution
3302    * @p __x into the output stream @p __os.
3303    *
3304    * @param __os An output stream.
3305    * @param __x  A %geometric_distribution random number distribution.
3306    *
3307    * @returns The output stream with the state of @p __x inserted or in
3308    * an error state.
3309    */
3310   template<typename _IntType,
3311            typename _CharT, typename _Traits>
3312     std::basic_ostream<_CharT, _Traits>&
3313     operator<<(std::basic_ostream<_CharT, _Traits>&,
3314                const std::geometric_distribution<_IntType>&);
3315
3316   /**
3317    * @brief Extracts a %geometric_distribution random number distribution
3318    * @p __x from the input stream @p __is.
3319    *
3320    * @param __is An input stream.
3321    * @param __x  A %geometric_distribution random number generator engine.
3322    *
3323    * @returns The input stream with @p __x extracted or in an error state.
3324    */
3325   template<typename _IntType,
3326            typename _CharT, typename _Traits>
3327     std::basic_istream<_CharT, _Traits>&
3328     operator>>(std::basic_istream<_CharT, _Traits>&,
3329                std::geometric_distribution<_IntType>&);
3330
3331
3332   /**
3333    * @brief A negative_binomial_distribution random number distribution.
3334    *
3335    * The formula for the negative binomial probability mass function is
3336    * @f$p(i) = \binom{n}{i} p^i (1 - p)^{t - i}@f$ where @f$t@f$
3337    * and @f$p@f$ are the parameters of the distribution.
3338    */
3339   template<typename _IntType = int>
3340     class negative_binomial_distribution
3341     {
3342       static_assert(std::is_integral<_IntType>::value,
3343                     "template argument not an integral type");
3344
3345     public:
3346       /** The type of the range of the distribution. */
3347       typedef _IntType result_type;
3348       /** Parameter type. */
3349       struct param_type
3350       {
3351         typedef negative_binomial_distribution<_IntType> distribution_type;
3352
3353         explicit
3354         param_type(_IntType __k = 1, double __p = 0.5)
3355         : _M_k(__k), _M_p(__p)
3356         { }
3357
3358         _IntType
3359         k() const
3360         { return _M_k; }
3361
3362         double
3363         p() const
3364         { return _M_p; }
3365
3366       private:
3367         _IntType _M_k;
3368         double _M_p;
3369       };
3370
3371       explicit
3372       negative_binomial_distribution(_IntType __k = 1, double __p = 0.5)
3373       : _M_param(__k, __p), _M_gd(__k, __p / (1.0 - __p))
3374       { }
3375
3376       explicit
3377       negative_binomial_distribution(const param_type& __p)
3378       : _M_param(__p), _M_gd(__p.k(), __p.p() / (1.0 - __p.p()))
3379       { }
3380
3381       /**
3382        * @brief Resets the distribution state.
3383        */
3384       void
3385       reset()
3386       { _M_gd.reset(); }
3387
3388       /**
3389        * @brief Return the @f$k@f$ parameter of the distribution.
3390        */
3391       _IntType
3392       k() const
3393       { return _M_param.k(); }
3394
3395       /**
3396        * @brief Return the @f$p@f$ parameter of the distribution.
3397        */
3398       double
3399       p() const
3400       { return _M_param.p(); }
3401
3402       /**
3403        * @brief Returns the parameter set of the distribution.
3404        */
3405       param_type
3406       param() const
3407       { return _M_param; }
3408
3409       /**
3410        * @brief Sets the parameter set of the distribution.
3411        * @param __param The new parameter set of the distribution.
3412        */
3413       void
3414       param(const param_type& __param)
3415       { _M_param = __param; }
3416
3417       /**
3418        * @brief Returns the greatest lower bound value of the distribution.
3419        */
3420       result_type
3421       min() const
3422       { return result_type(0); }
3423
3424       /**
3425        * @brief Returns the least upper bound value of the distribution.
3426        */
3427       result_type
3428       max() const
3429       { return std::numeric_limits<result_type>::max(); }
3430
3431       template<typename _UniformRandomNumberGenerator>
3432         result_type
3433         operator()(_UniformRandomNumberGenerator& __urng);
3434
3435       template<typename _UniformRandomNumberGenerator>
3436         result_type
3437         operator()(_UniformRandomNumberGenerator& __urng,
3438                    const param_type& __p);
3439
3440       /**
3441        * @brief Inserts a %negative_binomial_distribution random
3442        *        number distribution @p __x into the output stream @p __os.
3443        *
3444        * @param __os An output stream.
3445        * @param __x  A %negative_binomial_distribution random number
3446        *             distribution.
3447        *
3448        * @returns The output stream with the state of @p __x inserted or in
3449        *          an error state.
3450        */
3451       template<typename _IntType1, typename _CharT, typename _Traits>
3452         friend std::basic_ostream<_CharT, _Traits>&
3453         operator<<(std::basic_ostream<_CharT, _Traits>&,
3454                    const std::negative_binomial_distribution<_IntType1>&);
3455
3456       /**
3457        * @brief Extracts a %negative_binomial_distribution random number
3458        *        distribution @p __x from the input stream @p __is.
3459        *
3460        * @param __is An input stream.
3461        * @param __x A %negative_binomial_distribution random number
3462        *            generator engine.
3463        *
3464        * @returns The input stream with @p __x extracted or in an error state.
3465        */
3466       template<typename _IntType1, typename _CharT, typename _Traits>
3467         friend std::basic_istream<_CharT, _Traits>&
3468         operator>>(std::basic_istream<_CharT, _Traits>&,
3469                    std::negative_binomial_distribution<_IntType1>&);
3470
3471     private:
3472       param_type _M_param;
3473
3474       std::gamma_distribution<double> _M_gd;
3475     };
3476
3477   /* @} */ // group random_distributions_bernoulli
3478
3479   /**
3480    * @addtogroup random_distributions_poisson Poisson Distributions
3481    * @ingroup random_distributions
3482    * @{
3483    */
3484
3485   /**
3486    * @brief A discrete Poisson random number distribution.
3487    *
3488    * The formula for the Poisson probability density function is
3489    * @f$p(i|\mu) = \frac{\mu^i}{i!} e^{-\mu}@f$ where @f$\mu@f$ is the
3490    * parameter of the distribution.
3491    */
3492   template<typename _IntType = int>
3493     class poisson_distribution
3494     {
3495       static_assert(std::is_integral<_IntType>::value,
3496                     "template argument not an integral type");
3497
3498     public:
3499       /** The type of the range of the distribution. */
3500       typedef _IntType  result_type;
3501       /** Parameter type. */
3502       struct param_type
3503       {
3504         typedef poisson_distribution<_IntType> distribution_type;
3505         friend class poisson_distribution<_IntType>;
3506
3507         explicit
3508         param_type(double __mean = 1.0)
3509         : _M_mean(__mean)
3510         {
3511           _GLIBCXX_DEBUG_ASSERT(_M_mean > 0.0);
3512           _M_initialize();
3513         }
3514
3515         double
3516         mean() const
3517         { return _M_mean; }
3518
3519       private:
3520         // Hosts either log(mean) or the threshold of the simple method.
3521         void
3522         _M_initialize();
3523
3524         double _M_mean;
3525
3526         double _M_lm_thr;
3527 #if _GLIBCXX_USE_C99_MATH_TR1
3528         double _M_lfm, _M_sm, _M_d, _M_scx, _M_1cx, _M_c2b, _M_cb;
3529 #endif
3530       };
3531
3532       // constructors and member function
3533       explicit
3534       poisson_distribution(double __mean = 1.0)
3535       : _M_param(__mean), _M_nd()
3536       { }
3537
3538       explicit
3539       poisson_distribution(const param_type& __p)
3540       : _M_param(__p), _M_nd()
3541       { }
3542
3543       /**
3544        * @brief Resets the distribution state.
3545        */
3546       void
3547       reset()
3548       { _M_nd.reset(); }
3549
3550       /**
3551        * @brief Returns the distribution parameter @p mean.
3552        */
3553       double
3554       mean() const
3555       { return _M_param.mean(); }
3556
3557       /**
3558        * @brief Returns the parameter set of the distribution.
3559        */
3560       param_type
3561       param() const
3562       { return _M_param; }
3563
3564       /**
3565        * @brief Sets the parameter set of the distribution.
3566        * @param __param The new parameter set of the distribution.
3567        */
3568       void
3569       param(const param_type& __param)
3570       { _M_param = __param; }
3571
3572       /**
3573        * @brief Returns the greatest lower bound value of the distribution.
3574        */
3575       result_type
3576       min() const
3577       { return 0; }
3578
3579       /**
3580        * @brief Returns the least upper bound value of the distribution.
3581        */
3582       result_type
3583       max() const
3584       { return std::numeric_limits<result_type>::max(); }
3585
3586       template<typename _UniformRandomNumberGenerator>
3587         result_type
3588         operator()(_UniformRandomNumberGenerator& __urng)
3589         { return this->operator()(__urng, this->param()); }
3590
3591       template<typename _UniformRandomNumberGenerator>
3592         result_type
3593         operator()(_UniformRandomNumberGenerator& __urng,
3594                    const param_type& __p);
3595
3596       /**
3597        * @brief Inserts a %poisson_distribution random number distribution
3598        * @p __x into the output stream @p __os.
3599        *
3600        * @param __os An output stream.
3601        * @param __x  A %poisson_distribution random number distribution.
3602        *
3603        * @returns The output stream with the state of @p __x inserted or in
3604        * an error state.
3605        */
3606       template<typename _IntType1, typename _CharT, typename _Traits>
3607         friend std::basic_ostream<_CharT, _Traits>&
3608         operator<<(std::basic_ostream<_CharT, _Traits>&,
3609                    const std::poisson_distribution<_IntType1>&);
3610
3611       /**
3612        * @brief Extracts a %poisson_distribution random number distribution
3613        * @p __x from the input stream @p __is.
3614        *
3615        * @param __is An input stream.
3616        * @param __x  A %poisson_distribution random number generator engine.
3617        *
3618        * @returns The input stream with @p __x extracted or in an error
3619        *          state.
3620        */
3621       template<typename _IntType1, typename _CharT, typename _Traits>
3622         friend std::basic_istream<_CharT, _Traits>&
3623         operator>>(std::basic_istream<_CharT, _Traits>&,
3624                    std::poisson_distribution<_IntType1>&);
3625
3626     private:
3627       param_type _M_param;
3628
3629       // NB: Unused when _GLIBCXX_USE_C99_MATH_TR1 is undefined.
3630       std::normal_distribution<double> _M_nd;
3631     };
3632
3633   /**
3634    * @brief An exponential continuous distribution for random numbers.
3635    *
3636    * The formula for the exponential probability density function is
3637    * @f$p(x|\lambda) = \lambda e^{-\lambda x}@f$.
3638    *
3639    * <table border=1 cellpadding=10 cellspacing=0>
3640    * <caption align=top>Distribution Statistics</caption>
3641    * <tr><td>Mean</td><td>@f$\frac{1}{\lambda}@f$</td></tr>
3642    * <tr><td>Median</td><td>@f$\frac{\ln 2}{\lambda}@f$</td></tr>
3643    * <tr><td>Mode</td><td>@f$zero@f$</td></tr>
3644    * <tr><td>Range</td><td>@f$[0, \infty]@f$</td></tr>
3645    * <tr><td>Standard Deviation</td><td>@f$\frac{1}{\lambda}@f$</td></tr>
3646    * </table>
3647    */
3648   template<typename _RealType = double>
3649     class exponential_distribution
3650     {
3651       static_assert(std::is_floating_point<_RealType>::value,
3652                     "template argument not a floating point type");
3653
3654     public:
3655       /** The type of the range of the distribution. */
3656       typedef _RealType result_type;
3657       /** Parameter type. */
3658       struct param_type
3659       {
3660         typedef exponential_distribution<_RealType> distribution_type;
3661
3662         explicit
3663         param_type(_RealType __lambda = _RealType(1))
3664         : _M_lambda(__lambda)
3665         {
3666           _GLIBCXX_DEBUG_ASSERT(_M_lambda > _RealType(0));
3667         }
3668
3669         _RealType
3670         lambda() const
3671         { return _M_lambda; }
3672
3673       private:
3674         _RealType _M_lambda;
3675       };
3676
3677     public:
3678       /**
3679        * @brief Constructs an exponential distribution with inverse scale
3680        *        parameter @f$\lambda@f$.
3681        */
3682       explicit
3683       exponential_distribution(const result_type& __lambda = result_type(1))
3684       : _M_param(__lambda)
3685       { }
3686
3687       explicit
3688       exponential_distribution(const param_type& __p)
3689       : _M_param(__p)
3690       { }
3691
3692       /**
3693        * @brief Resets the distribution state.
3694        *
3695        * Has no effect on exponential distributions.
3696        */
3697       void
3698       reset() { }
3699
3700       /**
3701        * @brief Returns the inverse scale parameter of the distribution.
3702        */
3703       _RealType
3704       lambda() const
3705       { return _M_param.lambda(); }
3706
3707       /**
3708        * @brief Returns the parameter set of the distribution.
3709        */
3710       param_type
3711       param() const
3712       { return _M_param; }
3713
3714       /**
3715        * @brief Sets the parameter set of the distribution.
3716        * @param __param The new parameter set of the distribution.
3717        */
3718       void
3719       param(const param_type& __param)
3720       { _M_param = __param; }
3721
3722       /**
3723        * @brief Returns the greatest lower bound value of the distribution.
3724        */
3725       result_type
3726       min() const
3727       { return result_type(0); }
3728
3729       /**
3730        * @brief Returns the least upper bound value of the distribution.
3731        */
3732       result_type
3733       max() const
3734       { return std::numeric_limits<result_type>::max(); }
3735
3736       template<typename _UniformRandomNumberGenerator>
3737         result_type
3738         operator()(_UniformRandomNumberGenerator& __urng)
3739         { return this->operator()(__urng, this->param()); }
3740
3741       template<typename _UniformRandomNumberGenerator>
3742         result_type
3743         operator()(_UniformRandomNumberGenerator& __urng,
3744                    const param_type& __p)
3745         {
3746           __detail::_Adaptor<_UniformRandomNumberGenerator, result_type>
3747             __aurng(__urng);
3748           return -std::log(__aurng()) / __p.lambda();
3749         }
3750
3751     private:
3752       param_type _M_param;
3753     };
3754
3755   /**
3756    * @brief Inserts a %exponential_distribution random number distribution
3757    * @p __x into the output stream @p __os.
3758    *
3759    * @param __os An output stream.
3760    * @param __x  A %exponential_distribution random number distribution.
3761    *
3762    * @returns The output stream with the state of @p __x inserted or in
3763    * an error state.
3764    */
3765   template<typename _RealType, typename _CharT, typename _Traits>
3766     std::basic_ostream<_CharT, _Traits>&
3767     operator<<(std::basic_ostream<_CharT, _Traits>&,
3768                const std::exponential_distribution<_RealType>&);
3769
3770   /**
3771    * @brief Extracts a %exponential_distribution random number distribution
3772    * @p __x from the input stream @p __is.
3773    *
3774    * @param __is An input stream.
3775    * @param __x A %exponential_distribution random number
3776    *            generator engine.
3777    *
3778    * @returns The input stream with @p __x extracted or in an error state.
3779    */
3780   template<typename _RealType, typename _CharT, typename _Traits>
3781     std::basic_istream<_CharT, _Traits>&
3782     operator>>(std::basic_istream<_CharT, _Traits>&,
3783                std::exponential_distribution<_RealType>&);
3784
3785
3786   /**
3787    * @brief A weibull_distribution random number distribution.
3788    *
3789    * The formula for the normal probability density function is:
3790    * @f[
3791    *     p(x|\alpha,\beta) = \frac{\alpha}{\beta} (\frac{x}{\beta})^{\alpha-1}
3792    *                         \exp{(-(\frac{x}{\beta})^\alpha)}
3793    * @f]
3794    */
3795   template<typename _RealType = double>
3796     class weibull_distribution
3797     {
3798       static_assert(std::is_floating_point<_RealType>::value,
3799                     "template argument not a floating point type");
3800
3801     public:
3802       /** The type of the range of the distribution. */
3803       typedef _RealType result_type;
3804       /** Parameter type. */
3805       struct param_type
3806       {
3807         typedef weibull_distribution<_RealType> distribution_type;
3808
3809         explicit
3810         param_type(_RealType __a = _RealType(1),
3811                    _RealType __b = _RealType(1))
3812         : _M_a(__a), _M_b(__b)
3813         { }
3814
3815         _RealType
3816         a() const
3817         { return _M_a; }
3818
3819         _RealType
3820         b() const
3821         { return _M_b; }
3822
3823       private:
3824         _RealType _M_a;
3825         _RealType _M_b;
3826       };
3827
3828       explicit
3829       weibull_distribution(_RealType __a = _RealType(1),
3830                            _RealType __b = _RealType(1))
3831       : _M_param(__a, __b)
3832       { }
3833
3834       explicit
3835       weibull_distribution(const param_type& __p)
3836       : _M_param(__p)
3837       { }
3838
3839       /**
3840        * @brief Resets the distribution state.
3841        */
3842       void
3843       reset()
3844       { }
3845
3846       /**
3847        * @brief Return the @f$a@f$ parameter of the distribution.
3848        */
3849       _RealType
3850       a() const
3851       { return _M_param.a(); }
3852
3853       /**
3854        * @brief Return the @f$b@f$ parameter of the distribution.
3855        */
3856       _RealType
3857       b() const
3858       { return _M_param.b(); }
3859
3860       /**
3861        * @brief Returns the parameter set of the distribution.
3862        */
3863       param_type
3864       param() const
3865       { return _M_param; }
3866
3867       /**
3868        * @brief Sets the parameter set of the distribution.
3869        * @param __param The new parameter set of the distribution.
3870        */
3871       void
3872       param(const param_type& __param)
3873       { _M_param = __param; }
3874
3875       /**
3876        * @brief Returns the greatest lower bound value of the distribution.
3877        */
3878       result_type
3879       min() const
3880       { return result_type(0); }
3881
3882       /**
3883        * @brief Returns the least upper bound value of the distribution.
3884        */
3885       result_type
3886       max() const
3887       { return std::numeric_limits<result_type>::max(); }
3888
3889       template<typename _UniformRandomNumberGenerator>
3890         result_type
3891         operator()(_UniformRandomNumberGenerator& __urng)
3892         { return this->operator()(__urng, this->param()); }
3893
3894       template<typename _UniformRandomNumberGenerator>
3895         result_type
3896         operator()(_UniformRandomNumberGenerator& __urng,
3897                    const param_type& __p);
3898
3899     private:
3900       param_type _M_param;
3901     };
3902
3903   /**
3904    * @brief Inserts a %weibull_distribution random number distribution
3905    * @p __x into the output stream @p __os.
3906    *
3907    * @param __os An output stream.
3908    * @param __x  A %weibull_distribution random number distribution.
3909    *
3910    * @returns The output stream with the state of @p __x inserted or in
3911    * an error state.
3912    */
3913   template<typename _RealType, typename _CharT, typename _Traits>
3914     std::basic_ostream<_CharT, _Traits>&
3915     operator<<(std::basic_ostream<_CharT, _Traits>&,
3916                const std::weibull_distribution<_RealType>&);
3917
3918   /**
3919    * @brief Extracts a %weibull_distribution random number distribution
3920    * @p __x from the input stream @p __is.
3921    *
3922    * @param __is An input stream.
3923    * @param __x A %weibull_distribution random number
3924    *            generator engine.
3925    *
3926    * @returns The input stream with @p __x extracted or in an error state.
3927    */
3928   template<typename _RealType, typename _CharT, typename _Traits>
3929     std::basic_istream<_CharT, _Traits>&
3930     operator>>(std::basic_istream<_CharT, _Traits>&,
3931                std::weibull_distribution<_RealType>&);
3932
3933
3934   /**
3935    * @brief A extreme_value_distribution random number distribution.
3936    *
3937    * The formula for the normal probability mass function is
3938    * @f[
3939    *     p(x|a,b) = \frac{1}{b}
3940    *                \exp( \frac{a-x}{b} - \exp(\frac{a-x}{b}))
3941    * @f]
3942    */
3943   template<typename _RealType = double>
3944     class extreme_value_distribution
3945     {
3946       static_assert(std::is_floating_point<_RealType>::value,
3947                     "template argument not a floating point type");
3948
3949     public:
3950       /** The type of the range of the distribution. */
3951       typedef _RealType result_type;
3952       /** Parameter type. */
3953       struct param_type
3954       {
3955         typedef extreme_value_distribution<_RealType> distribution_type;
3956
3957         explicit
3958         param_type(_RealType __a = _RealType(0),
3959                    _RealType __b = _RealType(1))
3960         : _M_a(__a), _M_b(__b)
3961         { }
3962
3963         _RealType
3964         a() const
3965         { return _M_a; }
3966
3967         _RealType
3968         b() const
3969         { return _M_b; }
3970
3971       private:
3972         _RealType _M_a;
3973         _RealType _M_b;
3974       };
3975
3976       explicit
3977       extreme_value_distribution(_RealType __a = _RealType(0),
3978                                  _RealType __b = _RealType(1))
3979       : _M_param(__a, __b)
3980       { }
3981
3982       explicit
3983       extreme_value_distribution(const param_type& __p)
3984       : _M_param(__p)
3985       { }
3986
3987       /**
3988        * @brief Resets the distribution state.
3989        */
3990       void
3991       reset()
3992       { }
3993
3994       /**
3995        * @brief Return the @f$a@f$ parameter of the distribution.
3996        */
3997       _RealType
3998       a() const
3999       { return _M_param.a(); }
4000
4001       /**
4002        * @brief Return the @f$b@f$ parameter of the distribution.
4003        */
4004       _RealType
4005       b() const
4006       { return _M_param.b(); }
4007
4008       /**
4009        * @brief Returns the parameter set of the distribution.
4010        */
4011       param_type
4012       param() const
4013       { return _M_param; }
4014
4015       /**
4016        * @brief Sets the parameter set of the distribution.
4017        * @param __param The new parameter set of the distribution.
4018        */
4019       void
4020       param(const param_type& __param)
4021       { _M_param = __param; }
4022
4023       /**
4024        * @brief Returns the greatest lower bound value of the distribution.
4025        */
4026       result_type
4027       min() const
4028       { return std::numeric_limits<result_type>::min(); }
4029
4030       /**
4031        * @brief Returns the least upper bound value of the distribution.
4032        */
4033       result_type
4034       max() const
4035       { return std::numeric_limits<result_type>::max(); }
4036
4037       template<typename _UniformRandomNumberGenerator>
4038         result_type
4039         operator()(_UniformRandomNumberGenerator& __urng)
4040         { return this->operator()(__urng, this->param()); }
4041
4042       template<typename _UniformRandomNumberGenerator>
4043         result_type
4044         operator()(_UniformRandomNumberGenerator& __urng,
4045                    const param_type& __p);
4046
4047     private:
4048       param_type _M_param;
4049     };
4050
4051   /**
4052    * @brief Inserts a %extreme_value_distribution random number distribution
4053    * @p __x into the output stream @p __os.
4054    *
4055    * @param __os An output stream.
4056    * @param __x  A %extreme_value_distribution random number distribution.
4057    *
4058    * @returns The output stream with the state of @p __x inserted or in
4059    * an error state.
4060    */
4061   template<typename _RealType, typename _CharT, typename _Traits>
4062     std::basic_ostream<_CharT, _Traits>&
4063     operator<<(std::basic_ostream<_CharT, _Traits>&,
4064                const std::extreme_value_distribution<_RealType>&);
4065
4066   /**
4067    * @brief Extracts a %extreme_value_distribution random number
4068    *        distribution @p __x from the input stream @p __is.
4069    *
4070    * @param __is An input stream.
4071    * @param __x A %extreme_value_distribution random number
4072    *            generator engine.
4073    *
4074    * @returns The input stream with @p __x extracted or in an error state.
4075    */
4076   template<typename _RealType, typename _CharT, typename _Traits>
4077     std::basic_istream<_CharT, _Traits>&
4078     operator>>(std::basic_istream<_CharT, _Traits>&,
4079                std::extreme_value_distribution<_RealType>&);
4080
4081
4082   /**
4083    * @brief A discrete_distribution random number distribution.
4084    *
4085    * The formula for the discrete probability mass function is
4086    *
4087    */
4088   template<typename _IntType = int>
4089     class discrete_distribution
4090     {
4091       static_assert(std::is_integral<_IntType>::value,
4092                     "template argument not an integral type");
4093
4094     public:
4095       /** The type of the range of the distribution. */
4096       typedef _IntType result_type;
4097       /** Parameter type. */
4098       struct param_type
4099       {
4100         typedef discrete_distribution<_IntType> distribution_type;
4101         friend class discrete_distribution<_IntType>;
4102
4103         param_type()
4104         : _M_prob(), _M_cp()
4105         { _M_initialize(); }
4106
4107         template<typename _InputIterator>
4108           param_type(_InputIterator __wbegin,
4109                      _InputIterator __wend)
4110           : _M_prob(__wbegin, __wend), _M_cp()
4111           { _M_initialize(); }
4112
4113         param_type(initializer_list<double> __wil)
4114         : _M_prob(__wil.begin(), __wil.end()), _M_cp()
4115         { _M_initialize(); }
4116
4117         template<typename _Func>
4118           param_type(size_t __nw, double __xmin, double __xmax,
4119                      _Func __fw);
4120
4121         std::vector<double>
4122         probabilities() const
4123         { return _M_prob; }
4124
4125       private:
4126         void
4127         _M_initialize();
4128
4129         std::vector<double> _M_prob;
4130         std::vector<double> _M_cp;
4131       };
4132
4133       discrete_distribution()
4134       : _M_param()
4135       { }
4136
4137       template<typename _InputIterator>
4138         discrete_distribution(_InputIterator __wbegin,
4139                               _InputIterator __wend)
4140         : _M_param(__wbegin, __wend)
4141         { }
4142
4143       discrete_distribution(initializer_list<double> __wl)
4144       : _M_param(__wl)
4145       { }
4146
4147       template<typename _Func>
4148         discrete_distribution(size_t __nw, double __xmin, double __xmax,
4149                               _Func __fw)
4150         : _M_param(__nw, __xmin, __xmax, __fw)
4151         { }
4152
4153       explicit
4154       discrete_distribution(const param_type& __p)
4155       : _M_param(__p)
4156       { }
4157
4158       /**
4159        * @brief Resets the distribution state.
4160        */
4161       void
4162       reset()
4163       { }
4164
4165       /**
4166        * @brief Returns the probabilities of the distribution.
4167        */
4168       std::vector<double>
4169       probabilities() const
4170       { return _M_param.probabilities(); }
4171
4172       /**
4173        * @brief Returns the parameter set of the distribution.
4174        */
4175       param_type
4176       param() const
4177       { return _M_param; }
4178
4179       /**
4180        * @brief Sets the parameter set of the distribution.
4181        * @param __param The new parameter set of the distribution.
4182        */
4183       void
4184       param(const param_type& __param)
4185       { _M_param = __param; }
4186
4187       /**
4188        * @brief Returns the greatest lower bound value of the distribution.
4189        */
4190       result_type
4191       min() const
4192       { return result_type(0); }
4193
4194       /**
4195        * @brief Returns the least upper bound value of the distribution.
4196        */
4197       result_type
4198       max() const
4199       { return this->_M_param._M_prob.size() - 1; }
4200
4201       template<typename _UniformRandomNumberGenerator>
4202         result_type
4203         operator()(_UniformRandomNumberGenerator& __urng)
4204         { return this->operator()(__urng, this->param()); }
4205
4206       template<typename _UniformRandomNumberGenerator>
4207         result_type
4208         operator()(_UniformRandomNumberGenerator& __urng,
4209                    const param_type& __p);
4210
4211       /**
4212        * @brief Inserts a %discrete_distribution random number distribution
4213        * @p __x into the output stream @p __os.
4214        *
4215        * @param __os An output stream.
4216        * @param __x  A %discrete_distribution random number distribution.
4217        *
4218        * @returns The output stream with the state of @p __x inserted or in
4219        * an error state.
4220        */
4221       template<typename _IntType1, typename _CharT, typename _Traits>
4222         friend std::basic_ostream<_CharT, _Traits>&
4223         operator<<(std::basic_ostream<_CharT, _Traits>&,
4224                    const std::discrete_distribution<_IntType1>&);
4225
4226       /**
4227        * @brief Extracts a %discrete_distribution random number distribution
4228        * @p __x from the input stream @p __is.
4229        *
4230        * @param __is An input stream.
4231        * @param __x A %discrete_distribution random number
4232        *            generator engine.
4233        *
4234        * @returns The input stream with @p __x extracted or in an error
4235        *          state.
4236        */
4237       template<typename _IntType1, typename _CharT, typename _Traits>
4238         friend std::basic_istream<_CharT, _Traits>&
4239         operator>>(std::basic_istream<_CharT, _Traits>&,
4240                    std::discrete_distribution<_IntType1>&);
4241
4242     private:
4243       param_type _M_param;
4244     };
4245
4246
4247   /**
4248    * @brief A piecewise_constant_distribution random number distribution.
4249    *
4250    * The formula for the piecewise constant probability mass function is
4251    *
4252    */
4253   template<typename _RealType = double>
4254     class piecewise_constant_distribution
4255     {
4256       static_assert(std::is_floating_point<_RealType>::value,
4257                     "template argument not a floating point type");
4258
4259     public:
4260       /** The type of the range of the distribution. */
4261       typedef _RealType result_type;
4262       /** Parameter type. */
4263       struct param_type
4264       {
4265         typedef piecewise_constant_distribution<_RealType> distribution_type;
4266         friend class piecewise_constant_distribution<_RealType>;
4267
4268         param_type()
4269         : _M_int(), _M_den(), _M_cp()
4270         { _M_initialize(); }
4271
4272         template<typename _InputIteratorB, typename _InputIteratorW>
4273           param_type(_InputIteratorB __bfirst,
4274                      _InputIteratorB __bend,
4275                      _InputIteratorW __wbegin);
4276
4277         template<typename _Func>
4278           param_type(initializer_list<_RealType> __bi, _Func __fw);
4279
4280         template<typename _Func>
4281           param_type(size_t __nw, _RealType __xmin, _RealType __xmax,
4282                      _Func __fw);
4283
4284         std::vector<_RealType>
4285         intervals() const
4286         { return _M_int; }
4287
4288         std::vector<double>
4289         densities() const
4290         { return _M_den; }
4291
4292       private:
4293         void
4294         _M_initialize();
4295
4296         std::vector<_RealType> _M_int;
4297         std::vector<double> _M_den;
4298         std::vector<double> _M_cp;
4299       };
4300
4301       explicit
4302       piecewise_constant_distribution()
4303       : _M_param()
4304       { }
4305
4306       template<typename _InputIteratorB, typename _InputIteratorW>
4307         piecewise_constant_distribution(_InputIteratorB __bfirst,
4308                                         _InputIteratorB __bend,
4309                                         _InputIteratorW __wbegin)
4310         : _M_param(__bfirst, __bend, __wbegin)
4311         { }
4312
4313       template<typename _Func>
4314         piecewise_constant_distribution(initializer_list<_RealType> __bl,
4315                                         _Func __fw)
4316         : _M_param(__bl, __fw)
4317         { }
4318
4319       template<typename _Func>
4320         piecewise_constant_distribution(size_t __nw,
4321                                         _RealType __xmin, _RealType __xmax,
4322                                         _Func __fw)
4323         : _M_param(__nw, __xmin, __xmax, __fw)
4324         { }
4325
4326       explicit
4327       piecewise_constant_distribution(const param_type& __p)
4328       : _M_param(__p)
4329       { }
4330
4331       /**
4332        * @brief Resets the distribution state.
4333        */
4334       void
4335       reset()
4336       { }
4337
4338       /**
4339        * @brief Returns a vector of the intervals.
4340        */
4341       std::vector<_RealType>
4342       intervals() const
4343       { return _M_param.intervals(); }
4344
4345       /**
4346        * @brief Returns a vector of the probability densities.
4347        */
4348       std::vector<double>
4349       densities() const
4350       { return _M_param.densities(); }
4351
4352       /**
4353        * @brief Returns the parameter set of the distribution.
4354        */
4355       param_type
4356       param() const
4357       { return _M_param; }
4358
4359       /**
4360        * @brief Sets the parameter set of the distribution.
4361        * @param __param The new parameter set of the distribution.
4362        */
4363       void
4364       param(const param_type& __param)
4365       { _M_param = __param; }
4366
4367       /**
4368        * @brief Returns the greatest lower bound value of the distribution.
4369        */
4370       result_type
4371       min() const
4372       { return this->_M_param._M_int.front(); }
4373
4374       /**
4375        * @brief Returns the least upper bound value of the distribution.
4376        */
4377       result_type
4378       max() const
4379       { return this->_M_param._M_int.back(); }
4380
4381       template<typename _UniformRandomNumberGenerator>
4382         result_type
4383         operator()(_UniformRandomNumberGenerator& __urng)
4384         { return this->operator()(__urng, this->param()); }
4385
4386       template<typename _UniformRandomNumberGenerator>
4387         result_type
4388         operator()(_UniformRandomNumberGenerator& __urng,
4389                    const param_type& __p);
4390
4391       /**
4392        * @brief Inserts a %piecewise_constan_distribution random
4393        *        number distribution @p __x into the output stream @p __os.
4394        *
4395        * @param __os An output stream.
4396        * @param __x  A %piecewise_constan_distribution random number
4397        *             distribution.
4398        *
4399        * @returns The output stream with the state of @p __x inserted or in
4400        * an error state.
4401        */
4402       template<typename _RealType1, typename _CharT, typename _Traits>
4403         friend std::basic_ostream<_CharT, _Traits>&
4404         operator<<(std::basic_ostream<_CharT, _Traits>&,
4405                    const std::piecewise_constant_distribution<_RealType1>&);
4406
4407       /**
4408        * @brief Extracts a %piecewise_constan_distribution random
4409        *        number distribution @p __x from the input stream @p __is.
4410        *
4411        * @param __is An input stream.
4412        * @param __x A %piecewise_constan_distribution random number
4413        *            generator engine.
4414        *
4415        * @returns The input stream with @p __x extracted or in an error
4416        *          state.
4417        */
4418       template<typename _RealType1, typename _CharT, typename _Traits>
4419         friend std::basic_istream<_CharT, _Traits>&
4420         operator>>(std::basic_istream<_CharT, _Traits>&,
4421                    std::piecewise_constant_distribution<_RealType1>&);
4422
4423     private:
4424       param_type _M_param;
4425     };
4426
4427
4428   /**
4429    * @brief A piecewise_linear_distribution random number distribution.
4430    *
4431    * The formula for the piecewise linear probability mass function is
4432    *
4433    */
4434   template<typename _RealType = double>
4435     class piecewise_linear_distribution
4436     {
4437       static_assert(std::is_floating_point<_RealType>::value,
4438                     "template argument not a floating point type");
4439
4440     public:
4441       /** The type of the range of the distribution. */
4442       typedef _RealType result_type;
4443       /** Parameter type. */
4444       struct param_type
4445       {
4446         typedef piecewise_linear_distribution<_RealType> distribution_type;
4447         friend class piecewise_linear_distribution<_RealType>;
4448
4449         param_type()
4450         : _M_int(), _M_den(), _M_cp(), _M_m()
4451         { _M_initialize(); }
4452
4453         template<typename _InputIteratorB, typename _InputIteratorW>
4454           param_type(_InputIteratorB __bfirst,
4455                      _InputIteratorB __bend,
4456                      _InputIteratorW __wbegin);
4457
4458         template<typename _Func>
4459           param_type(initializer_list<_RealType> __bl, _Func __fw);
4460
4461         template<typename _Func>
4462           param_type(size_t __nw, _RealType __xmin, _RealType __xmax,
4463                      _Func __fw);
4464
4465         std::vector<_RealType>
4466         intervals() const
4467         { return _M_int; }
4468
4469         std::vector<double>
4470         densities() const
4471         { return _M_den; }
4472
4473       private:
4474         void
4475         _M_initialize();
4476
4477         std::vector<_RealType> _M_int;
4478         std::vector<double> _M_den;
4479         std::vector<double> _M_cp;
4480         std::vector<double> _M_m;
4481       };
4482
4483       explicit
4484       piecewise_linear_distribution()
4485       : _M_param()
4486       { }
4487
4488       template<typename _InputIteratorB, typename _InputIteratorW>
4489         piecewise_linear_distribution(_InputIteratorB __bfirst,
4490                                       _InputIteratorB __bend,
4491                                       _InputIteratorW __wbegin)
4492         : _M_param(__bfirst, __bend, __wbegin)
4493         { }
4494
4495       template<typename _Func>
4496         piecewise_linear_distribution(initializer_list<_RealType> __bl,
4497                                       _Func __fw)
4498         : _M_param(__bl, __fw)
4499         { }
4500
4501       template<typename _Func>
4502         piecewise_linear_distribution(size_t __nw,
4503                                       _RealType __xmin, _RealType __xmax,
4504                                       _Func __fw)
4505         : _M_param(__nw, __xmin, __xmax, __fw)
4506         { }
4507
4508       explicit
4509       piecewise_linear_distribution(const param_type& __p)
4510       : _M_param(__p)
4511       { }
4512
4513       /**
4514        * Resets the distribution state.
4515        */
4516       void
4517       reset()
4518       { }
4519
4520       /**
4521        * @brief Return the intervals of the distribution.
4522        */
4523       std::vector<_RealType>
4524       intervals() const
4525       { return _M_param.intervals(); }
4526
4527       /**
4528        * @brief Return a vector of the probability densities of the
4529        *        distribution.
4530        */
4531       std::vector<double>
4532       densities() const
4533       { return _M_param.densities(); }
4534
4535       /**
4536        * @brief Returns the parameter set of the distribution.
4537        */
4538       param_type
4539       param() const
4540       { return _M_param; }
4541
4542       /**
4543        * @brief Sets the parameter set of the distribution.
4544        * @param __param The new parameter set of the distribution.
4545        */
4546       void
4547       param(const param_type& __param)
4548       { _M_param = __param; }
4549
4550       /**
4551        * @brief Returns the greatest lower bound value of the distribution.
4552        */
4553       result_type
4554       min() const
4555       { return this->_M_param._M_int.front(); }
4556
4557       /**
4558        * @brief Returns the least upper bound value of the distribution.
4559        */
4560       result_type
4561       max() const
4562       { return this->_M_param._M_int.back(); }
4563
4564       template<typename _UniformRandomNumberGenerator>
4565         result_type
4566         operator()(_UniformRandomNumberGenerator& __urng)
4567         { return this->operator()(__urng, this->param()); }
4568
4569       template<typename _UniformRandomNumberGenerator>
4570         result_type
4571         operator()(_UniformRandomNumberGenerator& __urng,
4572                    const param_type& __p);
4573
4574       /**
4575        * @brief Inserts a %piecewise_linear_distribution random number
4576        *        distribution @p __x into the output stream @p __os.
4577        *
4578        * @param __os An output stream.
4579        * @param __x  A %piecewise_linear_distribution random number
4580        *             distribution.
4581        *
4582        * @returns The output stream with the state of @p __x inserted or in
4583        *          an error state.
4584        */
4585       template<typename _RealType1, typename _CharT, typename _Traits>
4586         friend std::basic_ostream<_CharT, _Traits>&
4587         operator<<(std::basic_ostream<_CharT, _Traits>&,
4588                    const std::piecewise_linear_distribution<_RealType1>&);
4589
4590       /**
4591        * @brief Extracts a %piecewise_linear_distribution random number
4592        *        distribution @p __x from the input stream @p __is.
4593        *
4594        * @param __is An input stream.
4595        * @param __x  A %piecewise_linear_distribution random number
4596        *             generator engine.
4597        *
4598        * @returns The input stream with @p __x extracted or in an error
4599        *          state.
4600        */
4601       template<typename _RealType1, typename _CharT, typename _Traits>
4602         friend std::basic_istream<_CharT, _Traits>&
4603         operator>>(std::basic_istream<_CharT, _Traits>&,
4604                    std::piecewise_linear_distribution<_RealType1>&);
4605
4606     private:
4607       param_type _M_param;
4608     };
4609
4610
4611   /* @} */ // group random_distributions_poisson
4612
4613   /* @} */ // group random_distributions
4614
4615   /**
4616    * @addtogroup random_utilities Random Number Utilities
4617    * @ingroup random
4618    * @{
4619    */
4620
4621   /**
4622    * @brief The seed_seq class generates sequences of seeds for random
4623    *        number generators.
4624    */
4625   class seed_seq
4626   {
4627
4628   public:
4629     /** The type of the seed vales. */
4630     typedef uint_least32_t result_type;
4631
4632     /** Default constructor. */
4633     seed_seq()
4634     : _M_v()
4635     { }
4636
4637     template<typename _IntType>
4638       seed_seq(std::initializer_list<_IntType> il);
4639
4640     template<typename _InputIterator>
4641       seed_seq(_InputIterator __begin, _InputIterator __end);
4642
4643     // generating functions
4644     template<typename _RandomAccessIterator>
4645       void
4646       generate(_RandomAccessIterator __begin, _RandomAccessIterator __end);
4647
4648     // property functions
4649     size_t size() const
4650     { return _M_v.size(); }
4651
4652     template<typename OutputIterator>
4653       void
4654       param(OutputIterator __dest) const
4655       { std::copy(_M_v.begin(), _M_v.end(), __dest); }
4656
4657   private:
4658     ///
4659     std::vector<result_type> _M_v;
4660   };
4661
4662   /* @} */ // group random_utilities
4663
4664   /* @} */ // group random
4665
4666 }
4667