Home | History | Annotate | Line # | Download | only in ext
random revision 1.1.1.3.2.1
      1 // Random number extensions -*- C++ -*-
      2 
      3 // Copyright (C) 2012-2017 Free Software Foundation, Inc.
      4 //
      5 // This file is part of the GNU ISO C++ Library.  This library is free
      6 // software; you can redistribute it and/or modify it under the
      7 // terms of the GNU General Public License as published by the
      8 // Free Software Foundation; either version 3, or (at your option)
      9 // any later version.
     10 
     11 // This library is distributed in the hope that it will be useful,
     12 // but WITHOUT ANY WARRANTY; without even the implied warranty of
     13 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
     14 // GNU General Public License for more details.
     15 
     16 // Under Section 7 of GPL version 3, you are granted additional
     17 // permissions described in the GCC Runtime Library Exception, version
     18 // 3.1, as published by the Free Software Foundation.
     19 
     20 // You should have received a copy of the GNU General Public License and
     21 // a copy of the GCC Runtime Library Exception along with this program;
     22 // see the files COPYING3 and COPYING.RUNTIME respectively.  If not, see
     23 // <http://www.gnu.org/licenses/>.
     24 
     25 /** @file ext/random
     26  *  This file is a GNU extension to the Standard C++ Library.
     27  */
     28 
     29 #ifndef _EXT_RANDOM
     30 #define _EXT_RANDOM 1
     31 
     32 #pragma GCC system_header
     33 
     34 #if __cplusplus < 201103L
     35 # include <bits/c++0x_warning.h>
     36 #else
     37 
     38 #include <random>
     39 #include <algorithm>
     40 #include <array>
     41 #include <ext/cmath>
     42 #ifdef __SSE2__
     43 # include <emmintrin.h>
     44 #endif
     45 
     46 #if defined(_GLIBCXX_USE_C99_STDINT_TR1) && defined(UINT32_C)
     47 
     48 namespace __gnu_cxx _GLIBCXX_VISIBILITY(default)
     49 {
     50 _GLIBCXX_BEGIN_NAMESPACE_VERSION
     51 
     52 #if __BYTE_ORDER__ == __ORDER_LITTLE_ENDIAN__
     53 
     54   /* Mersenne twister implementation optimized for vector operations.
     55    *
     56    * Reference: http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/SFMT/
     57    */
     58   template<typename _UIntType, size_t __m,
     59 	   size_t __pos1, size_t __sl1, size_t __sl2,
     60 	   size_t __sr1, size_t __sr2,
     61 	   uint32_t __msk1, uint32_t __msk2,
     62 	   uint32_t __msk3, uint32_t __msk4,
     63 	   uint32_t __parity1, uint32_t __parity2,
     64 	   uint32_t __parity3, uint32_t __parity4>
     65     class simd_fast_mersenne_twister_engine
     66     {
     67       static_assert(std::is_unsigned<_UIntType>::value, "template argument "
     68 		    "substituting _UIntType not an unsigned integral type");
     69       static_assert(__sr1 < 32, "first right shift too large");
     70       static_assert(__sr2 < 16, "second right shift too large");
     71       static_assert(__sl1 < 32, "first left shift too large");
     72       static_assert(__sl2 < 16, "second left shift too large");
     73 
     74     public:
     75       typedef _UIntType result_type;
     76 
     77     private:
     78       static constexpr size_t m_w = sizeof(result_type) * 8;
     79       static constexpr size_t _M_nstate = __m / 128 + 1;
     80       static constexpr size_t _M_nstate32 = _M_nstate * 4;
     81 
     82       static_assert(std::is_unsigned<_UIntType>::value, "template argument "
     83 		    "substituting _UIntType not an unsigned integral type");
     84       static_assert(__pos1 < _M_nstate, "POS1 not smaller than state size");
     85       static_assert(16 % sizeof(_UIntType) == 0,
     86 		    "UIntType size must divide 16");
     87 
     88     public:
     89       static constexpr size_t state_size = _M_nstate * (16
     90 							/ sizeof(result_type));
     91       static constexpr result_type default_seed = 5489u;
     92 
     93       // constructors and member function
     94       explicit
     95       simd_fast_mersenne_twister_engine(result_type __sd = default_seed)
     96       { seed(__sd); }
     97 
     98       template<typename _Sseq, typename = typename
     99 	std::enable_if<!std::is_same<_Sseq,
    100 				     simd_fast_mersenne_twister_engine>::value>
    101 	       ::type>
    102 	explicit
    103 	simd_fast_mersenne_twister_engine(_Sseq& __q)
    104 	{ seed(__q); }
    105 
    106       void
    107       seed(result_type __sd = default_seed);
    108 
    109       template<typename _Sseq>
    110 	typename std::enable_if<std::is_class<_Sseq>::value>::type
    111 	seed(_Sseq& __q);
    112 
    113       static constexpr result_type
    114       min()
    115       { return 0; };
    116 
    117       static constexpr result_type
    118       max()
    119       { return std::numeric_limits<result_type>::max(); }
    120 
    121       void
    122       discard(unsigned long long __z);
    123 
    124       result_type
    125       operator()()
    126       {
    127 	if (__builtin_expect(_M_pos >= state_size, 0))
    128 	  _M_gen_rand();
    129 
    130 	return _M_stateT[_M_pos++];
    131       }
    132 
    133       template<typename _UIntType_2, size_t __m_2,
    134 	       size_t __pos1_2, size_t __sl1_2, size_t __sl2_2,
    135 	       size_t __sr1_2, size_t __sr2_2,
    136 	       uint32_t __msk1_2, uint32_t __msk2_2,
    137 	       uint32_t __msk3_2, uint32_t __msk4_2,
    138 	       uint32_t __parity1_2, uint32_t __parity2_2,
    139 	       uint32_t __parity3_2, uint32_t __parity4_2>
    140 	friend bool
    141 	operator==(const simd_fast_mersenne_twister_engine<_UIntType_2,
    142 		   __m_2, __pos1_2, __sl1_2, __sl2_2, __sr1_2, __sr2_2,
    143 		   __msk1_2, __msk2_2, __msk3_2, __msk4_2,
    144 		   __parity1_2, __parity2_2, __parity3_2, __parity4_2>& __lhs,
    145 		   const simd_fast_mersenne_twister_engine<_UIntType_2,
    146 		   __m_2, __pos1_2, __sl1_2, __sl2_2, __sr1_2, __sr2_2,
    147 		   __msk1_2, __msk2_2, __msk3_2, __msk4_2,
    148 		   __parity1_2, __parity2_2, __parity3_2, __parity4_2>& __rhs);
    149 
    150       template<typename _UIntType_2, size_t __m_2,
    151 	       size_t __pos1_2, size_t __sl1_2, size_t __sl2_2,
    152 	       size_t __sr1_2, size_t __sr2_2,
    153 	       uint32_t __msk1_2, uint32_t __msk2_2,
    154 	       uint32_t __msk3_2, uint32_t __msk4_2,
    155 	       uint32_t __parity1_2, uint32_t __parity2_2,
    156 	       uint32_t __parity3_2, uint32_t __parity4_2,
    157 	       typename _CharT, typename _Traits>
    158 	friend std::basic_ostream<_CharT, _Traits>&
    159 	operator<<(std::basic_ostream<_CharT, _Traits>& __os,
    160 		   const __gnu_cxx::simd_fast_mersenne_twister_engine
    161 		   <_UIntType_2,
    162 		   __m_2, __pos1_2, __sl1_2, __sl2_2, __sr1_2, __sr2_2,
    163 		   __msk1_2, __msk2_2, __msk3_2, __msk4_2,
    164 		   __parity1_2, __parity2_2, __parity3_2, __parity4_2>& __x);
    165 
    166       template<typename _UIntType_2, size_t __m_2,
    167 	       size_t __pos1_2, size_t __sl1_2, size_t __sl2_2,
    168 	       size_t __sr1_2, size_t __sr2_2,
    169 	       uint32_t __msk1_2, uint32_t __msk2_2,
    170 	       uint32_t __msk3_2, uint32_t __msk4_2,
    171 	       uint32_t __parity1_2, uint32_t __parity2_2,
    172 	       uint32_t __parity3_2, uint32_t __parity4_2,
    173 	       typename _CharT, typename _Traits>
    174 	friend std::basic_istream<_CharT, _Traits>&
    175 	operator>>(std::basic_istream<_CharT, _Traits>& __is,
    176 		   __gnu_cxx::simd_fast_mersenne_twister_engine<_UIntType_2,
    177 		   __m_2, __pos1_2, __sl1_2, __sl2_2, __sr1_2, __sr2_2,
    178 		   __msk1_2, __msk2_2, __msk3_2, __msk4_2,
    179 		   __parity1_2, __parity2_2, __parity3_2, __parity4_2>& __x);
    180 
    181     private:
    182       union
    183       {
    184 #ifdef __SSE2__
    185 	__m128i _M_state[_M_nstate];
    186 #endif
    187 	uint32_t _M_state32[_M_nstate32];
    188 	result_type _M_stateT[state_size];
    189       } __attribute__ ((__aligned__ (16)));
    190       size_t _M_pos;
    191 
    192       void _M_gen_rand(void);
    193       void _M_period_certification();
    194   };
    195 
    196 
    197   template<typename _UIntType, size_t __m,
    198 	   size_t __pos1, size_t __sl1, size_t __sl2,
    199 	   size_t __sr1, size_t __sr2,
    200 	   uint32_t __msk1, uint32_t __msk2,
    201 	   uint32_t __msk3, uint32_t __msk4,
    202 	   uint32_t __parity1, uint32_t __parity2,
    203 	   uint32_t __parity3, uint32_t __parity4>
    204     inline bool
    205     operator!=(const __gnu_cxx::simd_fast_mersenne_twister_engine<_UIntType,
    206 	       __m, __pos1, __sl1, __sl2, __sr1, __sr2, __msk1, __msk2, __msk3,
    207 	       __msk4, __parity1, __parity2, __parity3, __parity4>& __lhs,
    208 	       const __gnu_cxx::simd_fast_mersenne_twister_engine<_UIntType,
    209 	       __m, __pos1, __sl1, __sl2, __sr1, __sr2, __msk1, __msk2, __msk3,
    210 	       __msk4, __parity1, __parity2, __parity3, __parity4>& __rhs)
    211     { return !(__lhs == __rhs); }
    212 
    213 
    214   /* Definitions for the SIMD-oriented Fast Mersenne Twister as defined
    215    * in the C implementation by Daito and Matsumoto, as both a 32-bit
    216    * and 64-bit version.
    217    */
    218   typedef simd_fast_mersenne_twister_engine<uint32_t, 607, 2,
    219 					    15, 3, 13, 3,
    220 					    0xfdff37ffU, 0xef7f3f7dU,
    221 					    0xff777b7dU, 0x7ff7fb2fU,
    222 					    0x00000001U, 0x00000000U,
    223 					    0x00000000U, 0x5986f054U>
    224     sfmt607;
    225 
    226   typedef simd_fast_mersenne_twister_engine<uint64_t, 607, 2,
    227 					    15, 3, 13, 3,
    228 					    0xfdff37ffU, 0xef7f3f7dU,
    229 					    0xff777b7dU, 0x7ff7fb2fU,
    230 					    0x00000001U, 0x00000000U,
    231 					    0x00000000U, 0x5986f054U>
    232     sfmt607_64;
    233 
    234 
    235   typedef simd_fast_mersenne_twister_engine<uint32_t, 1279, 7,
    236 					    14, 3, 5, 1,
    237 					    0xf7fefffdU, 0x7fefcfffU,
    238 					    0xaff3ef3fU, 0xb5ffff7fU,
    239 					    0x00000001U, 0x00000000U,
    240 					    0x00000000U, 0x20000000U>
    241     sfmt1279;
    242 
    243   typedef simd_fast_mersenne_twister_engine<uint64_t, 1279, 7,
    244 					    14, 3, 5, 1,
    245 					    0xf7fefffdU, 0x7fefcfffU,
    246 					    0xaff3ef3fU, 0xb5ffff7fU,
    247 					    0x00000001U, 0x00000000U,
    248 					    0x00000000U, 0x20000000U>
    249     sfmt1279_64;
    250 
    251 
    252   typedef simd_fast_mersenne_twister_engine<uint32_t, 2281, 12,
    253 					    19, 1, 5, 1,
    254 					    0xbff7ffbfU, 0xfdfffffeU,
    255 					    0xf7ffef7fU, 0xf2f7cbbfU,
    256 					    0x00000001U, 0x00000000U,
    257 					    0x00000000U, 0x41dfa600U>
    258     sfmt2281;
    259 
    260   typedef simd_fast_mersenne_twister_engine<uint64_t, 2281, 12,
    261 					    19, 1, 5, 1,
    262 					    0xbff7ffbfU, 0xfdfffffeU,
    263 					    0xf7ffef7fU, 0xf2f7cbbfU,
    264 					    0x00000001U, 0x00000000U,
    265 					    0x00000000U, 0x41dfa600U>
    266     sfmt2281_64;
    267 
    268 
    269   typedef simd_fast_mersenne_twister_engine<uint32_t, 4253, 17,
    270 					    20, 1, 7, 1,
    271 					    0x9f7bffffU, 0x9fffff5fU,
    272 					    0x3efffffbU, 0xfffff7bbU,
    273 					    0xa8000001U, 0xaf5390a3U,
    274 					    0xb740b3f8U, 0x6c11486dU>
    275     sfmt4253;
    276 
    277   typedef simd_fast_mersenne_twister_engine<uint64_t, 4253, 17,
    278 					    20, 1, 7, 1,
    279 					    0x9f7bffffU, 0x9fffff5fU,
    280 					    0x3efffffbU, 0xfffff7bbU,
    281 					    0xa8000001U, 0xaf5390a3U,
    282 					    0xb740b3f8U, 0x6c11486dU>
    283     sfmt4253_64;
    284 
    285 
    286   typedef simd_fast_mersenne_twister_engine<uint32_t, 11213, 68,
    287 					    14, 3, 7, 3,
    288 					    0xeffff7fbU, 0xffffffefU,
    289 					    0xdfdfbfffU, 0x7fffdbfdU,
    290 					    0x00000001U, 0x00000000U,
    291 					    0xe8148000U, 0xd0c7afa3U>
    292     sfmt11213;
    293 
    294   typedef simd_fast_mersenne_twister_engine<uint64_t, 11213, 68,
    295 					    14, 3, 7, 3,
    296 					    0xeffff7fbU, 0xffffffefU,
    297 					    0xdfdfbfffU, 0x7fffdbfdU,
    298 					    0x00000001U, 0x00000000U,
    299 					    0xe8148000U, 0xd0c7afa3U>
    300     sfmt11213_64;
    301 
    302 
    303   typedef simd_fast_mersenne_twister_engine<uint32_t, 19937, 122,
    304 					    18, 1, 11, 1,
    305 					    0xdfffffefU, 0xddfecb7fU,
    306 					    0xbffaffffU, 0xbffffff6U,
    307 					    0x00000001U, 0x00000000U,
    308 					    0x00000000U, 0x13c9e684U>
    309     sfmt19937;
    310 
    311   typedef simd_fast_mersenne_twister_engine<uint64_t, 19937, 122,
    312 					    18, 1, 11, 1,
    313 					    0xdfffffefU, 0xddfecb7fU,
    314 					    0xbffaffffU, 0xbffffff6U,
    315 					    0x00000001U, 0x00000000U,
    316 					    0x00000000U, 0x13c9e684U>
    317     sfmt19937_64;
    318 
    319 
    320   typedef simd_fast_mersenne_twister_engine<uint32_t, 44497, 330,
    321 					    5, 3, 9, 3,
    322 					    0xeffffffbU, 0xdfbebfffU,
    323 					    0xbfbf7befU, 0x9ffd7bffU,
    324 					    0x00000001U, 0x00000000U,
    325 					    0xa3ac4000U, 0xecc1327aU>
    326     sfmt44497;
    327 
    328   typedef simd_fast_mersenne_twister_engine<uint64_t, 44497, 330,
    329 					    5, 3, 9, 3,
    330 					    0xeffffffbU, 0xdfbebfffU,
    331 					    0xbfbf7befU, 0x9ffd7bffU,
    332 					    0x00000001U, 0x00000000U,
    333 					    0xa3ac4000U, 0xecc1327aU>
    334     sfmt44497_64;
    335 
    336 
    337   typedef simd_fast_mersenne_twister_engine<uint32_t, 86243, 366,
    338 					    6, 7, 19, 1,
    339 					    0xfdbffbffU, 0xbff7ff3fU,
    340 					    0xfd77efffU, 0xbf9ff3ffU,
    341 					    0x00000001U, 0x00000000U,
    342 					    0x00000000U, 0xe9528d85U>
    343     sfmt86243;
    344 
    345   typedef simd_fast_mersenne_twister_engine<uint64_t, 86243, 366,
    346 					    6, 7, 19, 1,
    347 					    0xfdbffbffU, 0xbff7ff3fU,
    348 					    0xfd77efffU, 0xbf9ff3ffU,
    349 					    0x00000001U, 0x00000000U,
    350 					    0x00000000U, 0xe9528d85U>
    351     sfmt86243_64;
    352 
    353 
    354   typedef simd_fast_mersenne_twister_engine<uint32_t, 132049, 110,
    355 					    19, 1, 21, 1,
    356 					    0xffffbb5fU, 0xfb6ebf95U,
    357 					    0xfffefffaU, 0xcff77fffU,
    358 					    0x00000001U, 0x00000000U,
    359 					    0xcb520000U, 0xc7e91c7dU>
    360     sfmt132049;
    361 
    362   typedef simd_fast_mersenne_twister_engine<uint64_t, 132049, 110,
    363 					    19, 1, 21, 1,
    364 					    0xffffbb5fU, 0xfb6ebf95U,
    365 					    0xfffefffaU, 0xcff77fffU,
    366 					    0x00000001U, 0x00000000U,
    367 					    0xcb520000U, 0xc7e91c7dU>
    368     sfmt132049_64;
    369 
    370 
    371   typedef simd_fast_mersenne_twister_engine<uint32_t, 216091, 627,
    372 					    11, 3, 10, 1,
    373 					    0xbff7bff7U, 0xbfffffffU,
    374 					    0xbffffa7fU, 0xffddfbfbU,
    375 					    0xf8000001U, 0x89e80709U,
    376 					    0x3bd2b64bU, 0x0c64b1e4U>
    377     sfmt216091;
    378 
    379   typedef simd_fast_mersenne_twister_engine<uint64_t, 216091, 627,
    380 					    11, 3, 10, 1,
    381 					    0xbff7bff7U, 0xbfffffffU,
    382 					    0xbffffa7fU, 0xffddfbfbU,
    383 					    0xf8000001U, 0x89e80709U,
    384 					    0x3bd2b64bU, 0x0c64b1e4U>
    385     sfmt216091_64;
    386 
    387 #endif // __BYTE_ORDER__ == __ORDER_LITTLE_ENDIAN__
    388 
    389   /**
    390    * @brief A beta continuous distribution for random numbers.
    391    *
    392    * The formula for the beta probability density function is:
    393    * @f[
    394    *     p(x|\alpha,\beta) = \frac{1}{B(\alpha,\beta)}
    395    *                         x^{\alpha - 1} (1 - x)^{\beta - 1}
    396    * @f]
    397    */
    398   template<typename _RealType = double>
    399     class beta_distribution
    400     {
    401       static_assert(std::is_floating_point<_RealType>::value,
    402 		    "template argument not a floating point type");
    403 
    404     public:
    405       /** The type of the range of the distribution. */
    406       typedef _RealType result_type;
    407 
    408       /** Parameter type. */
    409       struct param_type
    410       {
    411 	typedef beta_distribution<_RealType> distribution_type;
    412 	friend class beta_distribution<_RealType>;
    413 
    414 	explicit
    415 	param_type(_RealType __alpha_val = _RealType(1),
    416 		   _RealType __beta_val = _RealType(1))
    417 	: _M_alpha(__alpha_val), _M_beta(__beta_val)
    418 	{
    419 	  __glibcxx_assert(_M_alpha > _RealType(0));
    420 	  __glibcxx_assert(_M_beta > _RealType(0));
    421 	}
    422 
    423 	_RealType
    424 	alpha() const
    425 	{ return _M_alpha; }
    426 
    427 	_RealType
    428 	beta() const
    429 	{ return _M_beta; }
    430 
    431 	friend bool
    432 	operator==(const param_type& __p1, const param_type& __p2)
    433 	{ return (__p1._M_alpha == __p2._M_alpha
    434 		  && __p1._M_beta == __p2._M_beta); }
    435 
    436 	friend bool
    437 	operator!=(const param_type& __p1, const param_type& __p2)
    438 	{ return !(__p1 == __p2); }
    439 
    440       private:
    441 	void
    442 	_M_initialize();
    443 
    444 	_RealType _M_alpha;
    445 	_RealType _M_beta;
    446       };
    447 
    448     public:
    449       /**
    450        * @brief Constructs a beta distribution with parameters
    451        * @f$\alpha@f$ and @f$\beta@f$.
    452        */
    453       explicit
    454       beta_distribution(_RealType __alpha_val = _RealType(1),
    455 			_RealType __beta_val = _RealType(1))
    456       : _M_param(__alpha_val, __beta_val)
    457       { }
    458 
    459       explicit
    460       beta_distribution(const param_type& __p)
    461       : _M_param(__p)
    462       { }
    463 
    464       /**
    465        * @brief Resets the distribution state.
    466        */
    467       void
    468       reset()
    469       { }
    470 
    471       /**
    472        * @brief Returns the @f$\alpha@f$ of the distribution.
    473        */
    474       _RealType
    475       alpha() const
    476       { return _M_param.alpha(); }
    477 
    478       /**
    479        * @brief Returns the @f$\beta@f$ of the distribution.
    480        */
    481       _RealType
    482       beta() const
    483       { return _M_param.beta(); }
    484 
    485       /**
    486        * @brief Returns the parameter set of the distribution.
    487        */
    488       param_type
    489       param() const
    490       { return _M_param; }
    491 
    492       /**
    493        * @brief Sets the parameter set of the distribution.
    494        * @param __param The new parameter set of the distribution.
    495        */
    496       void
    497       param(const param_type& __param)
    498       { _M_param = __param; }
    499 
    500       /**
    501        * @brief Returns the greatest lower bound value of the distribution.
    502        */
    503       result_type
    504       min() const
    505       { return result_type(0); }
    506 
    507       /**
    508        * @brief Returns the least upper bound value of the distribution.
    509        */
    510       result_type
    511       max() const
    512       { return result_type(1); }
    513 
    514       /**
    515        * @brief Generating functions.
    516        */
    517       template<typename _UniformRandomNumberGenerator>
    518 	result_type
    519 	operator()(_UniformRandomNumberGenerator& __urng)
    520 	{ return this->operator()(__urng, _M_param); }
    521 
    522       template<typename _UniformRandomNumberGenerator>
    523 	result_type
    524 	operator()(_UniformRandomNumberGenerator& __urng,
    525 		   const param_type& __p);
    526 
    527       template<typename _ForwardIterator,
    528 	       typename _UniformRandomNumberGenerator>
    529 	void
    530 	__generate(_ForwardIterator __f, _ForwardIterator __t,
    531 		   _UniformRandomNumberGenerator& __urng)
    532 	{ this->__generate(__f, __t, __urng, _M_param); }
    533 
    534       template<typename _ForwardIterator,
    535 	       typename _UniformRandomNumberGenerator>
    536 	void
    537 	__generate(_ForwardIterator __f, _ForwardIterator __t,
    538 		   _UniformRandomNumberGenerator& __urng,
    539 		   const param_type& __p)
    540 	{ this->__generate_impl(__f, __t, __urng, __p); }
    541 
    542       template<typename _UniformRandomNumberGenerator>
    543 	void
    544 	__generate(result_type* __f, result_type* __t,
    545 		   _UniformRandomNumberGenerator& __urng,
    546 		   const param_type& __p)
    547 	{ this->__generate_impl(__f, __t, __urng, __p); }
    548 
    549       /**
    550        * @brief Return true if two beta distributions have the same
    551        *        parameters and the sequences that would be generated
    552        *        are equal.
    553        */
    554       friend bool
    555       operator==(const beta_distribution& __d1,
    556 		 const beta_distribution& __d2)
    557       { return __d1._M_param == __d2._M_param; }
    558 
    559       /**
    560        * @brief Inserts a %beta_distribution random number distribution
    561        * @p __x into the output stream @p __os.
    562        *
    563        * @param __os An output stream.
    564        * @param __x  A %beta_distribution random number distribution.
    565        *
    566        * @returns The output stream with the state of @p __x inserted or in
    567        * an error state.
    568        */
    569       template<typename _RealType1, typename _CharT, typename _Traits>
    570 	friend std::basic_ostream<_CharT, _Traits>&
    571 	operator<<(std::basic_ostream<_CharT, _Traits>& __os,
    572 		   const __gnu_cxx::beta_distribution<_RealType1>& __x);
    573 
    574       /**
    575        * @brief Extracts a %beta_distribution random number distribution
    576        * @p __x from the input stream @p __is.
    577        *
    578        * @param __is An input stream.
    579        * @param __x  A %beta_distribution random number generator engine.
    580        *
    581        * @returns The input stream with @p __x extracted or in an error state.
    582        */
    583       template<typename _RealType1, typename _CharT, typename _Traits>
    584 	friend std::basic_istream<_CharT, _Traits>&
    585 	operator>>(std::basic_istream<_CharT, _Traits>& __is,
    586 		   __gnu_cxx::beta_distribution<_RealType1>& __x);
    587 
    588     private:
    589       template<typename _ForwardIterator,
    590 	       typename _UniformRandomNumberGenerator>
    591 	void
    592 	__generate_impl(_ForwardIterator __f, _ForwardIterator __t,
    593 			_UniformRandomNumberGenerator& __urng,
    594 			const param_type& __p);
    595 
    596       param_type _M_param;
    597     };
    598 
    599   /**
    600    * @brief Return true if two beta distributions are different.
    601    */
    602   template<typename _RealType>
    603     inline bool
    604     operator!=(const __gnu_cxx::beta_distribution<_RealType>& __d1,
    605 	       const __gnu_cxx::beta_distribution<_RealType>& __d2)
    606     { return !(__d1 == __d2); }
    607 
    608 
    609   /**
    610    * @brief A multi-variate normal continuous distribution for random numbers.
    611    *
    612    * The formula for the normal probability density function is
    613    * @f[
    614    *     p(\overrightarrow{x}|\overrightarrow{\mu },\Sigma) =
    615    *       \frac{1}{\sqrt{(2\pi )^k\det(\Sigma))}}
    616    *       e^{-\frac{1}{2}(\overrightarrow{x}-\overrightarrow{\mu})^\text{T}
    617    *          \Sigma ^{-1}(\overrightarrow{x}-\overrightarrow{\mu})}
    618    * @f]
    619    *
    620    * where @f$\overrightarrow{x}@f$ and @f$\overrightarrow{\mu}@f$ are
    621    * vectors of dimension @f$k@f$ and @f$\Sigma@f$ is the covariance
    622    * matrix (which must be positive-definite).
    623    */
    624   template<std::size_t _Dimen, typename _RealType = double>
    625     class normal_mv_distribution
    626     {
    627       static_assert(std::is_floating_point<_RealType>::value,
    628 		    "template argument not a floating point type");
    629       static_assert(_Dimen != 0, "dimension is zero");
    630 
    631     public:
    632       /** The type of the range of the distribution. */
    633       typedef std::array<_RealType, _Dimen> result_type;
    634       /** Parameter type. */
    635       class param_type
    636       {
    637 	static constexpr size_t _M_t_size = _Dimen * (_Dimen + 1) / 2;
    638 
    639       public:
    640 	typedef normal_mv_distribution<_Dimen, _RealType> distribution_type;
    641 	friend class normal_mv_distribution<_Dimen, _RealType>;
    642 
    643 	param_type()
    644 	{
    645 	  std::fill(_M_mean.begin(), _M_mean.end(), _RealType(0));
    646 	  auto __it = _M_t.begin();
    647 	  for (size_t __i = 0; __i < _Dimen; ++__i)
    648 	    {
    649 	      std::fill_n(__it, __i, _RealType(0));
    650 	      __it += __i;
    651 	      *__it++ = _RealType(1);
    652 	    }
    653 	}
    654 
    655 	template<typename _ForwardIterator1, typename _ForwardIterator2>
    656 	  param_type(_ForwardIterator1 __meanbegin,
    657 		     _ForwardIterator1 __meanend,
    658 		     _ForwardIterator2 __varcovbegin,
    659 		     _ForwardIterator2 __varcovend)
    660 	{
    661 	  __glibcxx_function_requires(_ForwardIteratorConcept<
    662 				      _ForwardIterator1>)
    663 	  __glibcxx_function_requires(_ForwardIteratorConcept<
    664 				      _ForwardIterator2>)
    665 	  _GLIBCXX_DEBUG_ASSERT(std::distance(__meanbegin, __meanend)
    666 				<= _Dimen);
    667 	  const auto __dist = std::distance(__varcovbegin, __varcovend);
    668 	  _GLIBCXX_DEBUG_ASSERT(__dist == _Dimen * _Dimen
    669 				|| __dist == _Dimen * (_Dimen + 1) / 2
    670 				|| __dist == _Dimen);
    671 
    672 	  if (__dist == _Dimen * _Dimen)
    673 	    _M_init_full(__meanbegin, __meanend, __varcovbegin, __varcovend);
    674 	  else if (__dist == _Dimen * (_Dimen + 1) / 2)
    675 	    _M_init_lower(__meanbegin, __meanend, __varcovbegin, __varcovend);
    676 	  else
    677 	    {
    678 	      __glibcxx_assert(__dist == _Dimen);
    679 	      _M_init_diagonal(__meanbegin, __meanend,
    680 			       __varcovbegin, __varcovend);
    681 	    }
    682 	}
    683 
    684 	param_type(std::initializer_list<_RealType> __mean,
    685 		   std::initializer_list<_RealType> __varcov)
    686 	{
    687 	  _GLIBCXX_DEBUG_ASSERT(__mean.size() <= _Dimen);
    688 	  _GLIBCXX_DEBUG_ASSERT(__varcov.size() == _Dimen * _Dimen
    689 				|| __varcov.size() == _Dimen * (_Dimen + 1) / 2
    690 				|| __varcov.size() == _Dimen);
    691 
    692 	  if (__varcov.size() == _Dimen * _Dimen)
    693 	    _M_init_full(__mean.begin(), __mean.end(),
    694 			 __varcov.begin(), __varcov.end());
    695 	  else if (__varcov.size() == _Dimen * (_Dimen + 1) / 2)
    696 	    _M_init_lower(__mean.begin(), __mean.end(),
    697 			  __varcov.begin(), __varcov.end());
    698 	  else
    699 	    {
    700 	      __glibcxx_assert(__varcov.size() == _Dimen);
    701 	      _M_init_diagonal(__mean.begin(), __mean.end(),
    702 			       __varcov.begin(), __varcov.end());
    703 	    }
    704 	}
    705 
    706 	std::array<_RealType, _Dimen>
    707 	mean() const
    708 	{ return _M_mean; }
    709 
    710 	std::array<_RealType, _M_t_size>
    711 	varcov() const
    712 	{ return _M_t; }
    713 
    714 	friend bool
    715 	operator==(const param_type& __p1, const param_type& __p2)
    716 	{ return __p1._M_mean == __p2._M_mean && __p1._M_t == __p2._M_t; }
    717 
    718 	friend bool
    719 	operator!=(const param_type& __p1, const param_type& __p2)
    720 	{ return !(__p1 == __p2); }
    721 
    722       private:
    723 	template <typename _InputIterator1, typename _InputIterator2>
    724 	  void _M_init_full(_InputIterator1 __meanbegin,
    725 			    _InputIterator1 __meanend,
    726 			    _InputIterator2 __varcovbegin,
    727 			    _InputIterator2 __varcovend);
    728 	template <typename _InputIterator1, typename _InputIterator2>
    729 	  void _M_init_lower(_InputIterator1 __meanbegin,
    730 			     _InputIterator1 __meanend,
    731 			     _InputIterator2 __varcovbegin,
    732 			     _InputIterator2 __varcovend);
    733 	template <typename _InputIterator1, typename _InputIterator2>
    734 	  void _M_init_diagonal(_InputIterator1 __meanbegin,
    735 				_InputIterator1 __meanend,
    736 				_InputIterator2 __varbegin,
    737 				_InputIterator2 __varend);
    738 
    739 	std::array<_RealType, _Dimen> _M_mean;
    740 	std::array<_RealType, _M_t_size> _M_t;
    741       };
    742 
    743     public:
    744       normal_mv_distribution()
    745       : _M_param(), _M_nd()
    746       { }
    747 
    748       template<typename _ForwardIterator1, typename _ForwardIterator2>
    749 	normal_mv_distribution(_ForwardIterator1 __meanbegin,
    750 			       _ForwardIterator1 __meanend,
    751 			       _ForwardIterator2 __varcovbegin,
    752 			       _ForwardIterator2 __varcovend)
    753 	: _M_param(__meanbegin, __meanend, __varcovbegin, __varcovend),
    754 	  _M_nd()
    755 	{ }
    756 
    757       normal_mv_distribution(std::initializer_list<_RealType> __mean,
    758 			     std::initializer_list<_RealType> __varcov)
    759       : _M_param(__mean, __varcov), _M_nd()
    760       { }
    761 
    762       explicit
    763       normal_mv_distribution(const param_type& __p)
    764       : _M_param(__p), _M_nd()
    765       { }
    766 
    767       /**
    768        * @brief Resets the distribution state.
    769        */
    770       void
    771       reset()
    772       { _M_nd.reset(); }
    773 
    774       /**
    775        * @brief Returns the mean of the distribution.
    776        */
    777       result_type
    778       mean() const
    779       { return _M_param.mean(); }
    780 
    781       /**
    782        * @brief Returns the compact form of the variance/covariance
    783        * matrix of the distribution.
    784        */
    785       std::array<_RealType, _Dimen * (_Dimen + 1) / 2>
    786       varcov() const
    787       { return _M_param.varcov(); }
    788 
    789       /**
    790        * @brief Returns the parameter set of the distribution.
    791        */
    792       param_type
    793       param() const
    794       { return _M_param; }
    795 
    796       /**
    797        * @brief Sets the parameter set of the distribution.
    798        * @param __param The new parameter set of the distribution.
    799        */
    800       void
    801       param(const param_type& __param)
    802       { _M_param = __param; }
    803 
    804       /**
    805        * @brief Returns the greatest lower bound value of the distribution.
    806        */
    807       result_type
    808       min() const
    809       { result_type __res;
    810 	__res.fill(std::numeric_limits<_RealType>::lowest());
    811 	return __res; }
    812 
    813       /**
    814        * @brief Returns the least upper bound value of the distribution.
    815        */
    816       result_type
    817       max() const
    818       { result_type __res;
    819 	__res.fill(std::numeric_limits<_RealType>::max());
    820 	return __res; }
    821 
    822       /**
    823        * @brief Generating functions.
    824        */
    825       template<typename _UniformRandomNumberGenerator>
    826 	result_type
    827 	operator()(_UniformRandomNumberGenerator& __urng)
    828 	{ return this->operator()(__urng, _M_param); }
    829 
    830       template<typename _UniformRandomNumberGenerator>
    831 	result_type
    832 	operator()(_UniformRandomNumberGenerator& __urng,
    833 		   const param_type& __p);
    834 
    835       template<typename _ForwardIterator,
    836 	       typename _UniformRandomNumberGenerator>
    837 	void
    838 	__generate(_ForwardIterator __f, _ForwardIterator __t,
    839 		   _UniformRandomNumberGenerator& __urng)
    840 	{ return this->__generate_impl(__f, __t, __urng, _M_param); }
    841 
    842       template<typename _ForwardIterator,
    843 	       typename _UniformRandomNumberGenerator>
    844 	void
    845 	__generate(_ForwardIterator __f, _ForwardIterator __t,
    846 		   _UniformRandomNumberGenerator& __urng,
    847 		   const param_type& __p)
    848 	{ return this->__generate_impl(__f, __t, __urng, __p); }
    849 
    850       /**
    851        * @brief Return true if two multi-variant normal distributions have
    852        *        the same parameters and the sequences that would
    853        *        be generated are equal.
    854        */
    855       template<size_t _Dimen1, typename _RealType1>
    856 	friend bool
    857 	operator==(const
    858 		   __gnu_cxx::normal_mv_distribution<_Dimen1, _RealType1>&
    859 		   __d1,
    860 		   const
    861 		   __gnu_cxx::normal_mv_distribution<_Dimen1, _RealType1>&
    862 		   __d2);
    863 
    864       /**
    865        * @brief Inserts a %normal_mv_distribution random number distribution
    866        * @p __x into the output stream @p __os.
    867        *
    868        * @param __os An output stream.
    869        * @param __x  A %normal_mv_distribution random number distribution.
    870        *
    871        * @returns The output stream with the state of @p __x inserted or in
    872        * an error state.
    873        */
    874       template<size_t _Dimen1, typename _RealType1,
    875 	       typename _CharT, typename _Traits>
    876 	friend std::basic_ostream<_CharT, _Traits>&
    877 	operator<<(std::basic_ostream<_CharT, _Traits>& __os,
    878 		   const
    879 		   __gnu_cxx::normal_mv_distribution<_Dimen1, _RealType1>&
    880 		   __x);
    881 
    882       /**
    883        * @brief Extracts a %normal_mv_distribution random number distribution
    884        * @p __x from the input stream @p __is.
    885        *
    886        * @param __is An input stream.
    887        * @param __x  A %normal_mv_distribution random number generator engine.
    888        *
    889        * @returns The input stream with @p __x extracted or in an error
    890        *          state.
    891        */
    892       template<size_t _Dimen1, typename _RealType1,
    893 	       typename _CharT, typename _Traits>
    894 	friend std::basic_istream<_CharT, _Traits>&
    895 	operator>>(std::basic_istream<_CharT, _Traits>& __is,
    896 		   __gnu_cxx::normal_mv_distribution<_Dimen1, _RealType1>&
    897 		   __x);
    898 
    899     private:
    900       template<typename _ForwardIterator,
    901 	       typename _UniformRandomNumberGenerator>
    902 	void
    903 	__generate_impl(_ForwardIterator __f, _ForwardIterator __t,
    904 			_UniformRandomNumberGenerator& __urng,
    905 			const param_type& __p);
    906 
    907       param_type _M_param;
    908       std::normal_distribution<_RealType> _M_nd;
    909   };
    910 
    911   /**
    912    * @brief Return true if two multi-variate normal distributions are
    913    * different.
    914    */
    915   template<size_t _Dimen, typename _RealType>
    916     inline bool
    917     operator!=(const __gnu_cxx::normal_mv_distribution<_Dimen, _RealType>&
    918 	       __d1,
    919 	       const __gnu_cxx::normal_mv_distribution<_Dimen, _RealType>&
    920 	       __d2)
    921     { return !(__d1 == __d2); }
    922 
    923 
    924   /**
    925    * @brief A Rice continuous distribution for random numbers.
    926    *
    927    * The formula for the Rice probability density function is
    928    * @f[
    929    *     p(x|\nu,\sigma) = \frac{x}{\sigma^2}
    930    *                       \exp\left(-\frac{x^2+\nu^2}{2\sigma^2}\right)
    931    *                       I_0\left(\frac{x \nu}{\sigma^2}\right)
    932    * @f]
    933    * where @f$I_0(z)@f$ is the modified Bessel function of the first kind
    934    * of order 0 and @f$\nu >= 0@f$ and @f$\sigma > 0@f$.
    935    *
    936    * <table border=1 cellpadding=10 cellspacing=0>
    937    * <caption align=top>Distribution Statistics</caption>
    938    * <tr><td>Mean</td><td>@f$\sqrt{\pi/2}L_{1/2}(-\nu^2/2\sigma^2)@f$</td></tr>
    939    * <tr><td>Variance</td><td>@f$2\sigma^2 + \nu^2
    940    *                   + (\pi\sigma^2/2)L^2_{1/2}(-\nu^2/2\sigma^2)@f$</td></tr>
    941    * <tr><td>Range</td><td>@f$[0, \infty)@f$</td></tr>
    942    * </table>
    943    * where @f$L_{1/2}(x)@f$ is the Laguerre polynomial of order 1/2.
    944    */
    945   template<typename _RealType = double>
    946     class
    947     rice_distribution
    948     {
    949       static_assert(std::is_floating_point<_RealType>::value,
    950 		    "template argument not a floating point type");
    951     public:
    952       /** The type of the range of the distribution. */
    953       typedef _RealType result_type;
    954 
    955       /** Parameter type. */
    956       struct param_type
    957       {
    958 	typedef rice_distribution<result_type> distribution_type;
    959 
    960 	param_type(result_type __nu_val = result_type(0),
    961 		   result_type __sigma_val = result_type(1))
    962 	: _M_nu(__nu_val), _M_sigma(__sigma_val)
    963 	{
    964 	  __glibcxx_assert(_M_nu >= result_type(0));
    965 	  __glibcxx_assert(_M_sigma > result_type(0));
    966 	}
    967 
    968 	result_type
    969 	nu() const
    970 	{ return _M_nu; }
    971 
    972 	result_type
    973 	sigma() const
    974 	{ return _M_sigma; }
    975 
    976 	friend bool
    977 	operator==(const param_type& __p1, const param_type& __p2)
    978 	{ return __p1._M_nu == __p2._M_nu && __p1._M_sigma == __p2._M_sigma; }
    979 
    980 	friend bool
    981 	operator!=(const param_type& __p1, const param_type& __p2)
    982 	{ return !(__p1 == __p2); }
    983 
    984       private:
    985 	void _M_initialize();
    986 
    987 	result_type _M_nu;
    988 	result_type _M_sigma;
    989       };
    990 
    991       /**
    992        * @brief Constructors.
    993        */
    994       explicit
    995       rice_distribution(result_type __nu_val = result_type(0),
    996 			result_type __sigma_val = result_type(1))
    997       : _M_param(__nu_val, __sigma_val),
    998 	_M_ndx(__nu_val, __sigma_val),
    999 	_M_ndy(result_type(0), __sigma_val)
   1000       { }
   1001 
   1002       explicit
   1003       rice_distribution(const param_type& __p)
   1004       : _M_param(__p),
   1005 	_M_ndx(__p.nu(), __p.sigma()),
   1006 	_M_ndy(result_type(0), __p.sigma())
   1007       { }
   1008 
   1009       /**
   1010        * @brief Resets the distribution state.
   1011        */
   1012       void
   1013       reset()
   1014       {
   1015 	_M_ndx.reset();
   1016 	_M_ndy.reset();
   1017       }
   1018 
   1019       /**
   1020        * @brief Return the parameters of the distribution.
   1021        */
   1022       result_type
   1023       nu() const
   1024       { return _M_param.nu(); }
   1025 
   1026       result_type
   1027       sigma() const
   1028       { return _M_param.sigma(); }
   1029 
   1030       /**
   1031        * @brief Returns the parameter set of the distribution.
   1032        */
   1033       param_type
   1034       param() const
   1035       { return _M_param; }
   1036 
   1037       /**
   1038        * @brief Sets the parameter set of the distribution.
   1039        * @param __param The new parameter set of the distribution.
   1040        */
   1041       void
   1042       param(const param_type& __param)
   1043       { _M_param = __param; }
   1044 
   1045       /**
   1046        * @brief Returns the greatest lower bound value of the distribution.
   1047        */
   1048       result_type
   1049       min() const
   1050       { return result_type(0); }
   1051 
   1052       /**
   1053        * @brief Returns the least upper bound value of the distribution.
   1054        */
   1055       result_type
   1056       max() const
   1057       { return std::numeric_limits<result_type>::max(); }
   1058 
   1059       /**
   1060        * @brief Generating functions.
   1061        */
   1062       template<typename _UniformRandomNumberGenerator>
   1063 	result_type
   1064 	operator()(_UniformRandomNumberGenerator& __urng)
   1065 	{
   1066 	  result_type __x = this->_M_ndx(__urng);
   1067 	  result_type __y = this->_M_ndy(__urng);
   1068 #if _GLIBCXX_USE_C99_MATH_TR1
   1069 	  return std::hypot(__x, __y);
   1070 #else
   1071 	  return std::sqrt(__x * __x + __y * __y);
   1072 #endif
   1073 	}
   1074 
   1075       template<typename _UniformRandomNumberGenerator>
   1076 	result_type
   1077 	operator()(_UniformRandomNumberGenerator& __urng,
   1078 		   const param_type& __p)
   1079 	{
   1080 	  typename std::normal_distribution<result_type>::param_type
   1081 	    __px(__p.nu(), __p.sigma()), __py(result_type(0), __p.sigma());
   1082 	  result_type __x = this->_M_ndx(__px, __urng);
   1083 	  result_type __y = this->_M_ndy(__py, __urng);
   1084 #if _GLIBCXX_USE_C99_MATH_TR1
   1085 	  return std::hypot(__x, __y);
   1086 #else
   1087 	  return std::sqrt(__x * __x + __y * __y);
   1088 #endif
   1089 	}
   1090 
   1091       template<typename _ForwardIterator,
   1092 	       typename _UniformRandomNumberGenerator>
   1093 	void
   1094 	__generate(_ForwardIterator __f, _ForwardIterator __t,
   1095 		   _UniformRandomNumberGenerator& __urng)
   1096 	{ this->__generate(__f, __t, __urng, _M_param); }
   1097 
   1098       template<typename _ForwardIterator,
   1099 	       typename _UniformRandomNumberGenerator>
   1100 	void
   1101 	__generate(_ForwardIterator __f, _ForwardIterator __t,
   1102 		   _UniformRandomNumberGenerator& __urng,
   1103 		   const param_type& __p)
   1104 	{ this->__generate_impl(__f, __t, __urng, __p); }
   1105 
   1106       template<typename _UniformRandomNumberGenerator>
   1107 	void
   1108 	__generate(result_type* __f, result_type* __t,
   1109 		   _UniformRandomNumberGenerator& __urng,
   1110 		   const param_type& __p)
   1111 	{ this->__generate_impl(__f, __t, __urng, __p); }
   1112 
   1113       /**
   1114        * @brief Return true if two Rice distributions have
   1115        *        the same parameters and the sequences that would
   1116        *        be generated are equal.
   1117        */
   1118       friend bool
   1119       operator==(const rice_distribution& __d1,
   1120 		 const rice_distribution& __d2)
   1121       { return (__d1._M_param == __d2._M_param
   1122 		&& __d1._M_ndx == __d2._M_ndx
   1123 		&& __d1._M_ndy == __d2._M_ndy); }
   1124 
   1125       /**
   1126        * @brief Inserts a %rice_distribution random number distribution
   1127        * @p __x into the output stream @p __os.
   1128        *
   1129        * @param __os An output stream.
   1130        * @param __x  A %rice_distribution random number distribution.
   1131        *
   1132        * @returns The output stream with the state of @p __x inserted or in
   1133        * an error state.
   1134        */
   1135       template<typename _RealType1, typename _CharT, typename _Traits>
   1136 	friend std::basic_ostream<_CharT, _Traits>&
   1137 	operator<<(std::basic_ostream<_CharT, _Traits>&,
   1138 		   const rice_distribution<_RealType1>&);
   1139 
   1140       /**
   1141        * @brief Extracts a %rice_distribution random number distribution
   1142        * @p __x from the input stream @p __is.
   1143        *
   1144        * @param __is An input stream.
   1145        * @param __x A %rice_distribution random number
   1146        *            generator engine.
   1147        *
   1148        * @returns The input stream with @p __x extracted or in an error state.
   1149        */
   1150       template<typename _RealType1, typename _CharT, typename _Traits>
   1151 	friend std::basic_istream<_CharT, _Traits>&
   1152 	operator>>(std::basic_istream<_CharT, _Traits>&,
   1153 		   rice_distribution<_RealType1>&);
   1154 
   1155     private:
   1156       template<typename _ForwardIterator,
   1157 	       typename _UniformRandomNumberGenerator>
   1158 	void
   1159 	__generate_impl(_ForwardIterator __f, _ForwardIterator __t,
   1160 			_UniformRandomNumberGenerator& __urng,
   1161 			const param_type& __p);
   1162 
   1163       param_type _M_param;
   1164 
   1165       std::normal_distribution<result_type> _M_ndx;
   1166       std::normal_distribution<result_type> _M_ndy;
   1167     };
   1168 
   1169   /**
   1170    * @brief Return true if two Rice distributions are not equal.
   1171    */
   1172   template<typename _RealType1>
   1173     inline bool
   1174     operator!=(const rice_distribution<_RealType1>& __d1,
   1175 	       const rice_distribution<_RealType1>& __d2)
   1176     { return !(__d1 == __d2); }
   1177 
   1178 
   1179   /**
   1180    * @brief A Nakagami continuous distribution for random numbers.
   1181    *
   1182    * The formula for the Nakagami probability density function is
   1183    * @f[
   1184    *     p(x|\mu,\omega) = \frac{2\mu^\mu}{\Gamma(\mu)\omega^\mu}
   1185    *                       x^{2\mu-1}e^{-\mu x / \omega}
   1186    * @f]
   1187    * where @f$\Gamma(z)@f$ is the gamma function and @f$\mu >= 0.5@f$
   1188    * and @f$\omega > 0@f$.
   1189    */
   1190   template<typename _RealType = double>
   1191     class
   1192     nakagami_distribution
   1193     {
   1194       static_assert(std::is_floating_point<_RealType>::value,
   1195 		    "template argument not a floating point type");
   1196 
   1197     public:
   1198       /** The type of the range of the distribution. */
   1199       typedef _RealType result_type;
   1200 
   1201       /** Parameter type. */
   1202       struct param_type
   1203       {
   1204 	typedef nakagami_distribution<result_type> distribution_type;
   1205 
   1206 	param_type(result_type __mu_val = result_type(1),
   1207 		   result_type __omega_val = result_type(1))
   1208 	: _M_mu(__mu_val), _M_omega(__omega_val)
   1209 	{
   1210 	  __glibcxx_assert(_M_mu >= result_type(0.5L));
   1211 	  __glibcxx_assert(_M_omega > result_type(0));
   1212 	}
   1213 
   1214 	result_type
   1215 	mu() const
   1216 	{ return _M_mu; }
   1217 
   1218 	result_type
   1219 	omega() const
   1220 	{ return _M_omega; }
   1221 
   1222 	friend bool
   1223 	operator==(const param_type& __p1, const param_type& __p2)
   1224 	{ return __p1._M_mu == __p2._M_mu && __p1._M_omega == __p2._M_omega; }
   1225 
   1226 	friend bool
   1227 	operator!=(const param_type& __p1, const param_type& __p2)
   1228 	{ return !(__p1 == __p2); }
   1229 
   1230       private:
   1231 	void _M_initialize();
   1232 
   1233 	result_type _M_mu;
   1234 	result_type _M_omega;
   1235       };
   1236 
   1237       /**
   1238        * @brief Constructors.
   1239        */
   1240       explicit
   1241       nakagami_distribution(result_type __mu_val = result_type(1),
   1242 			    result_type __omega_val = result_type(1))
   1243       : _M_param(__mu_val, __omega_val),
   1244 	_M_gd(__mu_val, __omega_val / __mu_val)
   1245       { }
   1246 
   1247       explicit
   1248       nakagami_distribution(const param_type& __p)
   1249       : _M_param(__p),
   1250 	_M_gd(__p.mu(), __p.omega() / __p.mu())
   1251       { }
   1252 
   1253       /**
   1254        * @brief Resets the distribution state.
   1255        */
   1256       void
   1257       reset()
   1258       { _M_gd.reset(); }
   1259 
   1260       /**
   1261        * @brief Return the parameters of the distribution.
   1262        */
   1263       result_type
   1264       mu() const
   1265       { return _M_param.mu(); }
   1266 
   1267       result_type
   1268       omega() const
   1269       { return _M_param.omega(); }
   1270 
   1271       /**
   1272        * @brief Returns the parameter set of the distribution.
   1273        */
   1274       param_type
   1275       param() const
   1276       { return _M_param; }
   1277 
   1278       /**
   1279        * @brief Sets the parameter set of the distribution.
   1280        * @param __param The new parameter set of the distribution.
   1281        */
   1282       void
   1283       param(const param_type& __param)
   1284       { _M_param = __param; }
   1285 
   1286       /**
   1287        * @brief Returns the greatest lower bound value of the distribution.
   1288        */
   1289       result_type
   1290       min() const
   1291       { return result_type(0); }
   1292 
   1293       /**
   1294        * @brief Returns the least upper bound value of the distribution.
   1295        */
   1296       result_type
   1297       max() const
   1298       { return std::numeric_limits<result_type>::max(); }
   1299 
   1300       /**
   1301        * @brief Generating functions.
   1302        */
   1303       template<typename _UniformRandomNumberGenerator>
   1304 	result_type
   1305 	operator()(_UniformRandomNumberGenerator& __urng)
   1306 	{ return std::sqrt(this->_M_gd(__urng)); }
   1307 
   1308       template<typename _UniformRandomNumberGenerator>
   1309 	result_type
   1310 	operator()(_UniformRandomNumberGenerator& __urng,
   1311 		   const param_type& __p)
   1312 	{
   1313 	  typename std::gamma_distribution<result_type>::param_type
   1314 	    __pg(__p.mu(), __p.omega() / __p.mu());
   1315 	  return std::sqrt(this->_M_gd(__pg, __urng));
   1316 	}
   1317 
   1318       template<typename _ForwardIterator,
   1319 	       typename _UniformRandomNumberGenerator>
   1320 	void
   1321 	__generate(_ForwardIterator __f, _ForwardIterator __t,
   1322 		   _UniformRandomNumberGenerator& __urng)
   1323 	{ this->__generate(__f, __t, __urng, _M_param); }
   1324 
   1325       template<typename _ForwardIterator,
   1326 	       typename _UniformRandomNumberGenerator>
   1327 	void
   1328 	__generate(_ForwardIterator __f, _ForwardIterator __t,
   1329 		   _UniformRandomNumberGenerator& __urng,
   1330 		   const param_type& __p)
   1331 	{ this->__generate_impl(__f, __t, __urng, __p); }
   1332 
   1333       template<typename _UniformRandomNumberGenerator>
   1334 	void
   1335 	__generate(result_type* __f, result_type* __t,
   1336 		   _UniformRandomNumberGenerator& __urng,
   1337 		   const param_type& __p)
   1338 	{ this->__generate_impl(__f, __t, __urng, __p); }
   1339 
   1340       /**
   1341        * @brief Return true if two Nakagami distributions have
   1342        *        the same parameters and the sequences that would
   1343        *        be generated are equal.
   1344        */
   1345       friend bool
   1346       operator==(const nakagami_distribution& __d1,
   1347 		 const nakagami_distribution& __d2)
   1348       { return (__d1._M_param == __d2._M_param
   1349 		&& __d1._M_gd == __d2._M_gd); }
   1350 
   1351       /**
   1352        * @brief Inserts a %nakagami_distribution random number distribution
   1353        * @p __x into the output stream @p __os.
   1354        *
   1355        * @param __os An output stream.
   1356        * @param __x  A %nakagami_distribution random number distribution.
   1357        *
   1358        * @returns The output stream with the state of @p __x inserted or in
   1359        * an error state.
   1360        */
   1361       template<typename _RealType1, typename _CharT, typename _Traits>
   1362 	friend std::basic_ostream<_CharT, _Traits>&
   1363 	operator<<(std::basic_ostream<_CharT, _Traits>&,
   1364 		   const nakagami_distribution<_RealType1>&);
   1365 
   1366       /**
   1367        * @brief Extracts a %nakagami_distribution random number distribution
   1368        * @p __x from the input stream @p __is.
   1369        *
   1370        * @param __is An input stream.
   1371        * @param __x A %nakagami_distribution random number
   1372        *            generator engine.
   1373        *
   1374        * @returns The input stream with @p __x extracted or in an error state.
   1375        */
   1376       template<typename _RealType1, typename _CharT, typename _Traits>
   1377 	friend std::basic_istream<_CharT, _Traits>&
   1378 	operator>>(std::basic_istream<_CharT, _Traits>&,
   1379 		   nakagami_distribution<_RealType1>&);
   1380 
   1381     private:
   1382       template<typename _ForwardIterator,
   1383 	       typename _UniformRandomNumberGenerator>
   1384 	void
   1385 	__generate_impl(_ForwardIterator __f, _ForwardIterator __t,
   1386 			_UniformRandomNumberGenerator& __urng,
   1387 			const param_type& __p);
   1388 
   1389       param_type _M_param;
   1390 
   1391       std::gamma_distribution<result_type> _M_gd;
   1392     };
   1393 
   1394   /**
   1395    * @brief Return true if two Nakagami distributions are not equal.
   1396    */
   1397   template<typename _RealType>
   1398     inline bool
   1399     operator!=(const nakagami_distribution<_RealType>& __d1,
   1400 	       const nakagami_distribution<_RealType>& __d2)
   1401     { return !(__d1 == __d2); }
   1402 
   1403 
   1404   /**
   1405    * @brief A Pareto continuous distribution for random numbers.
   1406    *
   1407    * The formula for the Pareto cumulative probability function is
   1408    * @f[
   1409    *     P(x|\alpha,\mu) = 1 - \left(\frac{\mu}{x}\right)^\alpha
   1410    * @f]
   1411    * The formula for the Pareto probability density function is
   1412    * @f[
   1413    *     p(x|\alpha,\mu) = \frac{\alpha + 1}{\mu}
   1414    *                       \left(\frac{\mu}{x}\right)^{\alpha + 1}
   1415    * @f]
   1416    * where @f$x >= \mu@f$ and @f$\mu > 0@f$, @f$\alpha > 0@f$.
   1417    *
   1418    * <table border=1 cellpadding=10 cellspacing=0>
   1419    * <caption align=top>Distribution Statistics</caption>
   1420    * <tr><td>Mean</td><td>@f$\alpha \mu / (\alpha - 1)@f$
   1421    *              for @f$\alpha > 1@f$</td></tr>
   1422    * <tr><td>Variance</td><td>@f$\alpha \mu^2 / [(\alpha - 1)^2(\alpha - 2)]@f$
   1423    *              for @f$\alpha > 2@f$</td></tr>
   1424    * <tr><td>Range</td><td>@f$[\mu, \infty)@f$</td></tr>
   1425    * </table>
   1426    */
   1427   template<typename _RealType = double>
   1428     class
   1429     pareto_distribution
   1430     {
   1431       static_assert(std::is_floating_point<_RealType>::value,
   1432 		    "template argument not a floating point type");
   1433 
   1434     public:
   1435       /** The type of the range of the distribution. */
   1436       typedef _RealType result_type;
   1437 
   1438       /** Parameter type. */
   1439       struct param_type
   1440       {
   1441 	typedef pareto_distribution<result_type> distribution_type;
   1442 
   1443 	param_type(result_type __alpha_val = result_type(1),
   1444 		   result_type __mu_val = result_type(1))
   1445 	: _M_alpha(__alpha_val), _M_mu(__mu_val)
   1446 	{
   1447 	  __glibcxx_assert(_M_alpha > result_type(0));
   1448 	  __glibcxx_assert(_M_mu > result_type(0));
   1449 	}
   1450 
   1451 	result_type
   1452 	alpha() const
   1453 	{ return _M_alpha; }
   1454 
   1455 	result_type
   1456 	mu() const
   1457 	{ return _M_mu; }
   1458 
   1459 	friend bool
   1460 	operator==(const param_type& __p1, const param_type& __p2)
   1461 	{ return __p1._M_alpha == __p2._M_alpha && __p1._M_mu == __p2._M_mu; }
   1462 
   1463 	friend bool
   1464 	operator!=(const param_type& __p1, const param_type& __p2)
   1465 	{ return !(__p1 == __p2); }
   1466 
   1467       private:
   1468 	void _M_initialize();
   1469 
   1470 	result_type _M_alpha;
   1471 	result_type _M_mu;
   1472       };
   1473 
   1474       /**
   1475        * @brief Constructors.
   1476        */
   1477       explicit
   1478       pareto_distribution(result_type __alpha_val = result_type(1),
   1479 			  result_type __mu_val = result_type(1))
   1480       : _M_param(__alpha_val, __mu_val),
   1481 	_M_ud()
   1482       { }
   1483 
   1484       explicit
   1485       pareto_distribution(const param_type& __p)
   1486       : _M_param(__p),
   1487 	_M_ud()
   1488       { }
   1489 
   1490       /**
   1491        * @brief Resets the distribution state.
   1492        */
   1493       void
   1494       reset()
   1495       {
   1496 	_M_ud.reset();
   1497       }
   1498 
   1499       /**
   1500        * @brief Return the parameters of the distribution.
   1501        */
   1502       result_type
   1503       alpha() const
   1504       { return _M_param.alpha(); }
   1505 
   1506       result_type
   1507       mu() const
   1508       { return _M_param.mu(); }
   1509 
   1510       /**
   1511        * @brief Returns the parameter set of the distribution.
   1512        */
   1513       param_type
   1514       param() const
   1515       { return _M_param; }
   1516 
   1517       /**
   1518        * @brief Sets the parameter set of the distribution.
   1519        * @param __param The new parameter set of the distribution.
   1520        */
   1521       void
   1522       param(const param_type& __param)
   1523       { _M_param = __param; }
   1524 
   1525       /**
   1526        * @brief Returns the greatest lower bound value of the distribution.
   1527        */
   1528       result_type
   1529       min() const
   1530       { return this->mu(); }
   1531 
   1532       /**
   1533        * @brief Returns the least upper bound value of the distribution.
   1534        */
   1535       result_type
   1536       max() const
   1537       { return std::numeric_limits<result_type>::max(); }
   1538 
   1539       /**
   1540        * @brief Generating functions.
   1541        */
   1542       template<typename _UniformRandomNumberGenerator>
   1543 	result_type
   1544 	operator()(_UniformRandomNumberGenerator& __urng)
   1545 	{
   1546 	  return this->mu() * std::pow(this->_M_ud(__urng),
   1547 				       -result_type(1) / this->alpha());
   1548 	}
   1549 
   1550       template<typename _UniformRandomNumberGenerator>
   1551 	result_type
   1552 	operator()(_UniformRandomNumberGenerator& __urng,
   1553 		   const param_type& __p)
   1554 	{
   1555 	  return __p.mu() * std::pow(this->_M_ud(__urng),
   1556 					   -result_type(1) / __p.alpha());
   1557 	}
   1558 
   1559       template<typename _ForwardIterator,
   1560 	       typename _UniformRandomNumberGenerator>
   1561 	void
   1562 	__generate(_ForwardIterator __f, _ForwardIterator __t,
   1563 		   _UniformRandomNumberGenerator& __urng)
   1564 	{ this->__generate(__f, __t, __urng, _M_param); }
   1565 
   1566       template<typename _ForwardIterator,
   1567 	       typename _UniformRandomNumberGenerator>
   1568 	void
   1569 	__generate(_ForwardIterator __f, _ForwardIterator __t,
   1570 		   _UniformRandomNumberGenerator& __urng,
   1571 		   const param_type& __p)
   1572 	{ this->__generate_impl(__f, __t, __urng, __p); }
   1573 
   1574       template<typename _UniformRandomNumberGenerator>
   1575 	void
   1576 	__generate(result_type* __f, result_type* __t,
   1577 		   _UniformRandomNumberGenerator& __urng,
   1578 		   const param_type& __p)
   1579 	{ this->__generate_impl(__f, __t, __urng, __p); }
   1580 
   1581       /**
   1582        * @brief Return true if two Pareto distributions have
   1583        *        the same parameters and the sequences that would
   1584        *        be generated are equal.
   1585        */
   1586       friend bool
   1587       operator==(const pareto_distribution& __d1,
   1588 		 const pareto_distribution& __d2)
   1589       { return (__d1._M_param == __d2._M_param
   1590 		&& __d1._M_ud == __d2._M_ud); }
   1591 
   1592       /**
   1593        * @brief Inserts a %pareto_distribution random number distribution
   1594        * @p __x into the output stream @p __os.
   1595        *
   1596        * @param __os An output stream.
   1597        * @param __x  A %pareto_distribution random number distribution.
   1598        *
   1599        * @returns The output stream with the state of @p __x inserted or in
   1600        * an error state.
   1601        */
   1602       template<typename _RealType1, typename _CharT, typename _Traits>
   1603 	friend std::basic_ostream<_CharT, _Traits>&
   1604 	operator<<(std::basic_ostream<_CharT, _Traits>&,
   1605 		   const pareto_distribution<_RealType1>&);
   1606 
   1607       /**
   1608        * @brief Extracts a %pareto_distribution random number distribution
   1609        * @p __x from the input stream @p __is.
   1610        *
   1611        * @param __is An input stream.
   1612        * @param __x A %pareto_distribution random number
   1613        *            generator engine.
   1614        *
   1615        * @returns The input stream with @p __x extracted or in an error state.
   1616        */
   1617       template<typename _RealType1, typename _CharT, typename _Traits>
   1618 	friend std::basic_istream<_CharT, _Traits>&
   1619 	operator>>(std::basic_istream<_CharT, _Traits>&,
   1620 		   pareto_distribution<_RealType1>&);
   1621 
   1622     private:
   1623       template<typename _ForwardIterator,
   1624 	       typename _UniformRandomNumberGenerator>
   1625 	void
   1626 	__generate_impl(_ForwardIterator __f, _ForwardIterator __t,
   1627 			_UniformRandomNumberGenerator& __urng,
   1628 			const param_type& __p);
   1629 
   1630       param_type _M_param;
   1631 
   1632       std::uniform_real_distribution<result_type> _M_ud;
   1633     };
   1634 
   1635   /**
   1636    * @brief Return true if two Pareto distributions are not equal.
   1637    */
   1638   template<typename _RealType>
   1639     inline bool
   1640     operator!=(const pareto_distribution<_RealType>& __d1,
   1641 	       const pareto_distribution<_RealType>& __d2)
   1642     { return !(__d1 == __d2); }
   1643 
   1644 
   1645   /**
   1646    * @brief A K continuous distribution for random numbers.
   1647    *
   1648    * The formula for the K probability density function is
   1649    * @f[
   1650    *     p(x|\lambda, \mu, \nu) = \frac{2}{x}
   1651    *             \left(\frac{\lambda\nu x}{\mu}\right)^{\frac{\lambda + \nu}{2}}
   1652    *             \frac{1}{\Gamma(\lambda)\Gamma(\nu)}
   1653    *             K_{\nu - \lambda}\left(2\sqrt{\frac{\lambda\nu x}{\mu}}\right)
   1654    * @f]
   1655    * where @f$I_0(z)@f$ is the modified Bessel function of the second kind
   1656    * of order @f$\nu - \lambda@f$ and @f$\lambda > 0@f$, @f$\mu > 0@f$
   1657    * and @f$\nu > 0@f$.
   1658    *
   1659    * <table border=1 cellpadding=10 cellspacing=0>
   1660    * <caption align=top>Distribution Statistics</caption>
   1661    * <tr><td>Mean</td><td>@f$\mu@f$</td></tr>
   1662    * <tr><td>Variance</td><td>@f$\mu^2\frac{\lambda + \nu + 1}{\lambda\nu}@f$</td></tr>
   1663    * <tr><td>Range</td><td>@f$[0, \infty)@f$</td></tr>
   1664    * </table>
   1665    */
   1666   template<typename _RealType = double>
   1667     class
   1668     k_distribution
   1669     {
   1670       static_assert(std::is_floating_point<_RealType>::value,
   1671 		    "template argument not a floating point type");
   1672 
   1673     public:
   1674       /** The type of the range of the distribution. */
   1675       typedef _RealType result_type;
   1676 
   1677       /** Parameter type. */
   1678       struct param_type
   1679       {
   1680 	typedef k_distribution<result_type> distribution_type;
   1681 
   1682 	param_type(result_type __lambda_val = result_type(1),
   1683 		   result_type __mu_val = result_type(1),
   1684 		   result_type __nu_val = result_type(1))
   1685 	: _M_lambda(__lambda_val), _M_mu(__mu_val), _M_nu(__nu_val)
   1686 	{
   1687 	  __glibcxx_assert(_M_lambda > result_type(0));
   1688 	  __glibcxx_assert(_M_mu > result_type(0));
   1689 	  __glibcxx_assert(_M_nu > result_type(0));
   1690 	}
   1691 
   1692 	result_type
   1693 	lambda() const
   1694 	{ return _M_lambda; }
   1695 
   1696 	result_type
   1697 	mu() const
   1698 	{ return _M_mu; }
   1699 
   1700 	result_type
   1701 	nu() const
   1702 	{ return _M_nu; }
   1703 
   1704 	friend bool
   1705 	operator==(const param_type& __p1, const param_type& __p2)
   1706 	{
   1707 	  return __p1._M_lambda == __p2._M_lambda
   1708 	      && __p1._M_mu == __p2._M_mu
   1709 	      && __p1._M_nu == __p2._M_nu;
   1710 	}
   1711 
   1712 	friend bool
   1713 	operator!=(const param_type& __p1, const param_type& __p2)
   1714 	{ return !(__p1 == __p2); }
   1715 
   1716       private:
   1717 	void _M_initialize();
   1718 
   1719 	result_type _M_lambda;
   1720 	result_type _M_mu;
   1721 	result_type _M_nu;
   1722       };
   1723 
   1724       /**
   1725        * @brief Constructors.
   1726        */
   1727       explicit
   1728       k_distribution(result_type __lambda_val = result_type(1),
   1729 		     result_type __mu_val = result_type(1),
   1730 		     result_type __nu_val = result_type(1))
   1731       : _M_param(__lambda_val, __mu_val, __nu_val),
   1732 	_M_gd1(__lambda_val, result_type(1) / __lambda_val),
   1733 	_M_gd2(__nu_val, __mu_val / __nu_val)
   1734       { }
   1735 
   1736       explicit
   1737       k_distribution(const param_type& __p)
   1738       : _M_param(__p),
   1739 	_M_gd1(__p.lambda(), result_type(1) / __p.lambda()),
   1740 	_M_gd2(__p.nu(), __p.mu() / __p.nu())
   1741       { }
   1742 
   1743       /**
   1744        * @brief Resets the distribution state.
   1745        */
   1746       void
   1747       reset()
   1748       {
   1749 	_M_gd1.reset();
   1750 	_M_gd2.reset();
   1751       }
   1752 
   1753       /**
   1754        * @brief Return the parameters of the distribution.
   1755        */
   1756       result_type
   1757       lambda() const
   1758       { return _M_param.lambda(); }
   1759 
   1760       result_type
   1761       mu() const
   1762       { return _M_param.mu(); }
   1763 
   1764       result_type
   1765       nu() const
   1766       { return _M_param.nu(); }
   1767 
   1768       /**
   1769        * @brief Returns the parameter set of the distribution.
   1770        */
   1771       param_type
   1772       param() const
   1773       { return _M_param; }
   1774 
   1775       /**
   1776        * @brief Sets the parameter set of the distribution.
   1777        * @param __param The new parameter set of the distribution.
   1778        */
   1779       void
   1780       param(const param_type& __param)
   1781       { _M_param = __param; }
   1782 
   1783       /**
   1784        * @brief Returns the greatest lower bound value of the distribution.
   1785        */
   1786       result_type
   1787       min() const
   1788       { return result_type(0); }
   1789 
   1790       /**
   1791        * @brief Returns the least upper bound value of the distribution.
   1792        */
   1793       result_type
   1794       max() const
   1795       { return std::numeric_limits<result_type>::max(); }
   1796 
   1797       /**
   1798        * @brief Generating functions.
   1799        */
   1800       template<typename _UniformRandomNumberGenerator>
   1801 	result_type
   1802 	operator()(_UniformRandomNumberGenerator&);
   1803 
   1804       template<typename _UniformRandomNumberGenerator>
   1805 	result_type
   1806 	operator()(_UniformRandomNumberGenerator&, const param_type&);
   1807 
   1808       template<typename _ForwardIterator,
   1809 	       typename _UniformRandomNumberGenerator>
   1810 	void
   1811 	__generate(_ForwardIterator __f, _ForwardIterator __t,
   1812 		   _UniformRandomNumberGenerator& __urng)
   1813 	{ this->__generate(__f, __t, __urng, _M_param); }
   1814 
   1815       template<typename _ForwardIterator,
   1816 	       typename _UniformRandomNumberGenerator>
   1817 	void
   1818 	__generate(_ForwardIterator __f, _ForwardIterator __t,
   1819 		   _UniformRandomNumberGenerator& __urng,
   1820 		   const param_type& __p)
   1821 	{ this->__generate_impl(__f, __t, __urng, __p); }
   1822 
   1823       template<typename _UniformRandomNumberGenerator>
   1824 	void
   1825 	__generate(result_type* __f, result_type* __t,
   1826 		   _UniformRandomNumberGenerator& __urng,
   1827 		   const param_type& __p)
   1828 	{ this->__generate_impl(__f, __t, __urng, __p); }
   1829 
   1830       /**
   1831        * @brief Return true if two K distributions have
   1832        *        the same parameters and the sequences that would
   1833        *        be generated are equal.
   1834        */
   1835       friend bool
   1836       operator==(const k_distribution& __d1,
   1837 		 const k_distribution& __d2)
   1838       { return (__d1._M_param == __d2._M_param
   1839 		&& __d1._M_gd1 == __d2._M_gd1
   1840 		&& __d1._M_gd2 == __d2._M_gd2); }
   1841 
   1842       /**
   1843        * @brief Inserts a %k_distribution random number distribution
   1844        * @p __x into the output stream @p __os.
   1845        *
   1846        * @param __os An output stream.
   1847        * @param __x  A %k_distribution random number distribution.
   1848        *
   1849        * @returns The output stream with the state of @p __x inserted or in
   1850        * an error state.
   1851        */
   1852       template<typename _RealType1, typename _CharT, typename _Traits>
   1853 	friend std::basic_ostream<_CharT, _Traits>&
   1854 	operator<<(std::basic_ostream<_CharT, _Traits>&,
   1855 		   const k_distribution<_RealType1>&);
   1856 
   1857       /**
   1858        * @brief Extracts a %k_distribution random number distribution
   1859        * @p __x from the input stream @p __is.
   1860        *
   1861        * @param __is An input stream.
   1862        * @param __x A %k_distribution random number
   1863        *            generator engine.
   1864        *
   1865        * @returns The input stream with @p __x extracted or in an error state.
   1866        */
   1867       template<typename _RealType1, typename _CharT, typename _Traits>
   1868 	friend std::basic_istream<_CharT, _Traits>&
   1869 	operator>>(std::basic_istream<_CharT, _Traits>&,
   1870 		   k_distribution<_RealType1>&);
   1871 
   1872     private:
   1873       template<typename _ForwardIterator,
   1874 	       typename _UniformRandomNumberGenerator>
   1875 	void
   1876 	__generate_impl(_ForwardIterator __f, _ForwardIterator __t,
   1877 			_UniformRandomNumberGenerator& __urng,
   1878 			const param_type& __p);
   1879 
   1880       param_type _M_param;
   1881 
   1882       std::gamma_distribution<result_type> _M_gd1;
   1883       std::gamma_distribution<result_type> _M_gd2;
   1884     };
   1885 
   1886   /**
   1887    * @brief Return true if two K distributions are not equal.
   1888    */
   1889   template<typename _RealType>
   1890     inline bool
   1891     operator!=(const k_distribution<_RealType>& __d1,
   1892 	       const k_distribution<_RealType>& __d2)
   1893     { return !(__d1 == __d2); }
   1894 
   1895 
   1896   /**
   1897    * @brief An arcsine continuous distribution for random numbers.
   1898    *
   1899    * The formula for the arcsine probability density function is
   1900    * @f[
   1901    *     p(x|a,b) = \frac{1}{\pi \sqrt{(x - a)(b - x)}}
   1902    * @f]
   1903    * where @f$x >= a@f$ and @f$x <= b@f$.
   1904    *
   1905    * <table border=1 cellpadding=10 cellspacing=0>
   1906    * <caption align=top>Distribution Statistics</caption>
   1907    * <tr><td>Mean</td><td>@f$ (a + b) / 2 @f$</td></tr>
   1908    * <tr><td>Variance</td><td>@f$ (b - a)^2 / 8 @f$</td></tr>
   1909    * <tr><td>Range</td><td>@f$[a, b]@f$</td></tr>
   1910    * </table>
   1911    */
   1912   template<typename _RealType = double>
   1913     class
   1914     arcsine_distribution
   1915     {
   1916       static_assert(std::is_floating_point<_RealType>::value,
   1917 		    "template argument not a floating point type");
   1918 
   1919     public:
   1920       /** The type of the range of the distribution. */
   1921       typedef _RealType result_type;
   1922 
   1923       /** Parameter type. */
   1924       struct param_type
   1925       {
   1926 	typedef arcsine_distribution<result_type> distribution_type;
   1927 
   1928 	param_type(result_type __a = result_type(0),
   1929 		   result_type __b = result_type(1))
   1930 	: _M_a(__a), _M_b(__b)
   1931 	{
   1932 	  __glibcxx_assert(_M_a <= _M_b);
   1933 	}
   1934 
   1935 	result_type
   1936 	a() const
   1937 	{ return _M_a; }
   1938 
   1939 	result_type
   1940 	b() const
   1941 	{ return _M_b; }
   1942 
   1943 	friend bool
   1944 	operator==(const param_type& __p1, const param_type& __p2)
   1945 	{ return __p1._M_a == __p2._M_a && __p1._M_b == __p2._M_b; }
   1946 
   1947 	friend bool
   1948 	operator!=(const param_type& __p1, const param_type& __p2)
   1949 	{ return !(__p1 == __p2); }
   1950 
   1951       private:
   1952 	void _M_initialize();
   1953 
   1954 	result_type _M_a;
   1955 	result_type _M_b;
   1956       };
   1957 
   1958       /**
   1959        * @brief Constructors.
   1960        */
   1961       explicit
   1962       arcsine_distribution(result_type __a = result_type(0),
   1963 			   result_type __b = result_type(1))
   1964       : _M_param(__a, __b),
   1965 	_M_ud(-1.5707963267948966192313216916397514L,
   1966 	      +1.5707963267948966192313216916397514L)
   1967       { }
   1968 
   1969       explicit
   1970       arcsine_distribution(const param_type& __p)
   1971       : _M_param(__p),
   1972 	_M_ud(-1.5707963267948966192313216916397514L,
   1973 	      +1.5707963267948966192313216916397514L)
   1974       { }
   1975 
   1976       /**
   1977        * @brief Resets the distribution state.
   1978        */
   1979       void
   1980       reset()
   1981       { _M_ud.reset(); }
   1982 
   1983       /**
   1984        * @brief Return the parameters of the distribution.
   1985        */
   1986       result_type
   1987       a() const
   1988       { return _M_param.a(); }
   1989 
   1990       result_type
   1991       b() const
   1992       { return _M_param.b(); }
   1993 
   1994       /**
   1995        * @brief Returns the parameter set of the distribution.
   1996        */
   1997       param_type
   1998       param() const
   1999       { return _M_param; }
   2000 
   2001       /**
   2002        * @brief Sets the parameter set of the distribution.
   2003        * @param __param The new parameter set of the distribution.
   2004        */
   2005       void
   2006       param(const param_type& __param)
   2007       { _M_param = __param; }
   2008 
   2009       /**
   2010        * @brief Returns the greatest lower bound value of the distribution.
   2011        */
   2012       result_type
   2013       min() const
   2014       { return this->a(); }
   2015 
   2016       /**
   2017        * @brief Returns the least upper bound value of the distribution.
   2018        */
   2019       result_type
   2020       max() const
   2021       { return this->b(); }
   2022 
   2023       /**
   2024        * @brief Generating functions.
   2025        */
   2026       template<typename _UniformRandomNumberGenerator>
   2027 	result_type
   2028 	operator()(_UniformRandomNumberGenerator& __urng)
   2029 	{
   2030 	  result_type __x = std::sin(this->_M_ud(__urng));
   2031 	  return (__x * (this->b() - this->a())
   2032 		  + this->a() + this->b()) / result_type(2);
   2033 	}
   2034 
   2035       template<typename _UniformRandomNumberGenerator>
   2036 	result_type
   2037 	operator()(_UniformRandomNumberGenerator& __urng,
   2038 		   const param_type& __p)
   2039 	{
   2040 	  result_type __x = std::sin(this->_M_ud(__urng));
   2041 	  return (__x * (__p.b() - __p.a())
   2042 		  + __p.a() + __p.b()) / result_type(2);
   2043 	}
   2044 
   2045       template<typename _ForwardIterator,
   2046 	       typename _UniformRandomNumberGenerator>
   2047 	void
   2048 	__generate(_ForwardIterator __f, _ForwardIterator __t,
   2049 		   _UniformRandomNumberGenerator& __urng)
   2050 	{ this->__generate(__f, __t, __urng, _M_param); }
   2051 
   2052       template<typename _ForwardIterator,
   2053 	       typename _UniformRandomNumberGenerator>
   2054 	void
   2055 	__generate(_ForwardIterator __f, _ForwardIterator __t,
   2056 		   _UniformRandomNumberGenerator& __urng,
   2057 		   const param_type& __p)
   2058 	{ this->__generate_impl(__f, __t, __urng, __p); }
   2059 
   2060       template<typename _UniformRandomNumberGenerator>
   2061 	void
   2062 	__generate(result_type* __f, result_type* __t,
   2063 		   _UniformRandomNumberGenerator& __urng,
   2064 		   const param_type& __p)
   2065 	{ this->__generate_impl(__f, __t, __urng, __p); }
   2066 
   2067       /**
   2068        * @brief Return true if two arcsine distributions have
   2069        *        the same parameters and the sequences that would
   2070        *        be generated are equal.
   2071        */
   2072       friend bool
   2073       operator==(const arcsine_distribution& __d1,
   2074 		 const arcsine_distribution& __d2)
   2075       { return (__d1._M_param == __d2._M_param
   2076 		&& __d1._M_ud == __d2._M_ud); }
   2077 
   2078       /**
   2079        * @brief Inserts a %arcsine_distribution random number distribution
   2080        * @p __x into the output stream @p __os.
   2081        *
   2082        * @param __os An output stream.
   2083        * @param __x  A %arcsine_distribution random number distribution.
   2084        *
   2085        * @returns The output stream with the state of @p __x inserted or in
   2086        * an error state.
   2087        */
   2088       template<typename _RealType1, typename _CharT, typename _Traits>
   2089 	friend std::basic_ostream<_CharT, _Traits>&
   2090 	operator<<(std::basic_ostream<_CharT, _Traits>&,
   2091 		   const arcsine_distribution<_RealType1>&);
   2092 
   2093       /**
   2094        * @brief Extracts a %arcsine_distribution random number distribution
   2095        * @p __x from the input stream @p __is.
   2096        *
   2097        * @param __is An input stream.
   2098        * @param __x A %arcsine_distribution random number
   2099        *            generator engine.
   2100        *
   2101        * @returns The input stream with @p __x extracted or in an error state.
   2102        */
   2103       template<typename _RealType1, typename _CharT, typename _Traits>
   2104 	friend std::basic_istream<_CharT, _Traits>&
   2105 	operator>>(std::basic_istream<_CharT, _Traits>&,
   2106 		   arcsine_distribution<_RealType1>&);
   2107 
   2108     private:
   2109       template<typename _ForwardIterator,
   2110 	       typename _UniformRandomNumberGenerator>
   2111 	void
   2112 	__generate_impl(_ForwardIterator __f, _ForwardIterator __t,
   2113 			_UniformRandomNumberGenerator& __urng,
   2114 			const param_type& __p);
   2115 
   2116       param_type _M_param;
   2117 
   2118       std::uniform_real_distribution<result_type> _M_ud;
   2119     };
   2120 
   2121   /**
   2122    * @brief Return true if two arcsine distributions are not equal.
   2123    */
   2124   template<typename _RealType>
   2125     inline bool
   2126     operator!=(const arcsine_distribution<_RealType>& __d1,
   2127 	       const arcsine_distribution<_RealType>& __d2)
   2128     { return !(__d1 == __d2); }
   2129 
   2130 
   2131   /**
   2132    * @brief A Hoyt continuous distribution for random numbers.
   2133    *
   2134    * The formula for the Hoyt probability density function is
   2135    * @f[
   2136    *     p(x|q,\omega) = \frac{(1 + q^2)x}{q\omega}
   2137    *                     \exp\left(-\frac{(1 + q^2)^2 x^2}{4 q^2 \omega}\right)
   2138    *                       I_0\left(\frac{(1 - q^4) x^2}{4 q^2 \omega}\right)
   2139    * @f]
   2140    * where @f$I_0(z)@f$ is the modified Bessel function of the first kind
   2141    * of order 0 and @f$0 < q < 1@f$.
   2142    *
   2143    * <table border=1 cellpadding=10 cellspacing=0>
   2144    * <caption align=top>Distribution Statistics</caption>
   2145    * <tr><td>Mean</td><td>@f$ \sqrt{\frac{2}{\pi}} \sqrt{\frac{\omega}{1 + q^2}}
   2146    *                       E(1 - q^2) @f$</td></tr>
   2147    * <tr><td>Variance</td><td>@f$ \omega \left(1 - \frac{2E^2(1 - q^2)}
   2148    *                                      {\pi (1 + q^2)}\right) @f$</td></tr>
   2149    * <tr><td>Range</td><td>@f$[0, \infty)@f$</td></tr>
   2150    * </table>
   2151    * where @f$E(x)@f$ is the elliptic function of the second kind.
   2152    */
   2153   template<typename _RealType = double>
   2154     class
   2155     hoyt_distribution
   2156     {
   2157       static_assert(std::is_floating_point<_RealType>::value,
   2158 		    "template argument not a floating point type");
   2159 
   2160     public:
   2161       /** The type of the range of the distribution. */
   2162       typedef _RealType result_type;
   2163 
   2164       /** Parameter type. */
   2165       struct param_type
   2166       {
   2167 	typedef hoyt_distribution<result_type> distribution_type;
   2168 
   2169 	param_type(result_type __q = result_type(0.5L),
   2170 		   result_type __omega = result_type(1))
   2171 	: _M_q(__q), _M_omega(__omega)
   2172 	{
   2173 	  __glibcxx_assert(_M_q > result_type(0));
   2174 	  __glibcxx_assert(_M_q < result_type(1));
   2175 	}
   2176 
   2177 	result_type
   2178 	q() const
   2179 	{ return _M_q; }
   2180 
   2181 	result_type
   2182 	omega() const
   2183 	{ return _M_omega; }
   2184 
   2185 	friend bool
   2186 	operator==(const param_type& __p1, const param_type& __p2)
   2187 	{ return __p1._M_q == __p2._M_q && __p1._M_omega == __p2._M_omega; }
   2188 
   2189 	friend bool
   2190 	operator!=(const param_type& __p1, const param_type& __p2)
   2191 	{ return !(__p1 == __p2); }
   2192 
   2193       private:
   2194 	void _M_initialize();
   2195 
   2196 	result_type _M_q;
   2197 	result_type _M_omega;
   2198       };
   2199 
   2200       /**
   2201        * @brief Constructors.
   2202        */
   2203       explicit
   2204       hoyt_distribution(result_type __q = result_type(0.5L),
   2205 			result_type __omega = result_type(1))
   2206       : _M_param(__q, __omega),
   2207 	_M_ad(result_type(0.5L) * (result_type(1) + __q * __q),
   2208 	      result_type(0.5L) * (result_type(1) + __q * __q)
   2209 				/ (__q * __q)),
   2210 	_M_ed(result_type(1))
   2211       { }
   2212 
   2213       explicit
   2214       hoyt_distribution(const param_type& __p)
   2215       : _M_param(__p),
   2216 	_M_ad(result_type(0.5L) * (result_type(1) + __p.q() * __p.q()),
   2217 	      result_type(0.5L) * (result_type(1) + __p.q() * __p.q())
   2218 				/ (__p.q() * __p.q())),
   2219 	_M_ed(result_type(1))
   2220       { }
   2221 
   2222       /**
   2223        * @brief Resets the distribution state.
   2224        */
   2225       void
   2226       reset()
   2227       {
   2228 	_M_ad.reset();
   2229 	_M_ed.reset();
   2230       }
   2231 
   2232       /**
   2233        * @brief Return the parameters of the distribution.
   2234        */
   2235       result_type
   2236       q() const
   2237       { return _M_param.q(); }
   2238 
   2239       result_type
   2240       omega() const
   2241       { return _M_param.omega(); }
   2242 
   2243       /**
   2244        * @brief Returns the parameter set of the distribution.
   2245        */
   2246       param_type
   2247       param() const
   2248       { return _M_param; }
   2249 
   2250       /**
   2251        * @brief Sets the parameter set of the distribution.
   2252        * @param __param The new parameter set of the distribution.
   2253        */
   2254       void
   2255       param(const param_type& __param)
   2256       { _M_param = __param; }
   2257 
   2258       /**
   2259        * @brief Returns the greatest lower bound value of the distribution.
   2260        */
   2261       result_type
   2262       min() const
   2263       { return result_type(0); }
   2264 
   2265       /**
   2266        * @brief Returns the least upper bound value of the distribution.
   2267        */
   2268       result_type
   2269       max() const
   2270       { return std::numeric_limits<result_type>::max(); }
   2271 
   2272       /**
   2273        * @brief Generating functions.
   2274        */
   2275       template<typename _UniformRandomNumberGenerator>
   2276 	result_type
   2277 	operator()(_UniformRandomNumberGenerator& __urng);
   2278 
   2279       template<typename _UniformRandomNumberGenerator>
   2280 	result_type
   2281 	operator()(_UniformRandomNumberGenerator& __urng,
   2282 		   const param_type& __p);
   2283 
   2284       template<typename _ForwardIterator,
   2285 	       typename _UniformRandomNumberGenerator>
   2286 	void
   2287 	__generate(_ForwardIterator __f, _ForwardIterator __t,
   2288 		   _UniformRandomNumberGenerator& __urng)
   2289 	{ this->__generate(__f, __t, __urng, _M_param); }
   2290 
   2291       template<typename _ForwardIterator,
   2292 	       typename _UniformRandomNumberGenerator>
   2293 	void
   2294 	__generate(_ForwardIterator __f, _ForwardIterator __t,
   2295 		   _UniformRandomNumberGenerator& __urng,
   2296 		   const param_type& __p)
   2297 	{ this->__generate_impl(__f, __t, __urng, __p); }
   2298 
   2299       template<typename _UniformRandomNumberGenerator>
   2300 	void
   2301 	__generate(result_type* __f, result_type* __t,
   2302 		   _UniformRandomNumberGenerator& __urng,
   2303 		   const param_type& __p)
   2304 	{ this->__generate_impl(__f, __t, __urng, __p); }
   2305 
   2306       /**
   2307        * @brief Return true if two Hoyt distributions have
   2308        *        the same parameters and the sequences that would
   2309        *        be generated are equal.
   2310        */
   2311       friend bool
   2312       operator==(const hoyt_distribution& __d1,
   2313 		 const hoyt_distribution& __d2)
   2314       { return (__d1._M_param == __d2._M_param
   2315 		&& __d1._M_ad == __d2._M_ad
   2316 		&& __d1._M_ed == __d2._M_ed); }
   2317 
   2318       /**
   2319        * @brief Inserts a %hoyt_distribution random number distribution
   2320        * @p __x into the output stream @p __os.
   2321        *
   2322        * @param __os An output stream.
   2323        * @param __x  A %hoyt_distribution random number distribution.
   2324        *
   2325        * @returns The output stream with the state of @p __x inserted or in
   2326        * an error state.
   2327        */
   2328       template<typename _RealType1, typename _CharT, typename _Traits>
   2329 	friend std::basic_ostream<_CharT, _Traits>&
   2330 	operator<<(std::basic_ostream<_CharT, _Traits>&,
   2331 		   const hoyt_distribution<_RealType1>&);
   2332 
   2333       /**
   2334        * @brief Extracts a %hoyt_distribution random number distribution
   2335        * @p __x from the input stream @p __is.
   2336        *
   2337        * @param __is An input stream.
   2338        * @param __x A %hoyt_distribution random number
   2339        *            generator engine.
   2340        *
   2341        * @returns The input stream with @p __x extracted or in an error state.
   2342        */
   2343       template<typename _RealType1, typename _CharT, typename _Traits>
   2344 	friend std::basic_istream<_CharT, _Traits>&
   2345 	operator>>(std::basic_istream<_CharT, _Traits>&,
   2346 		   hoyt_distribution<_RealType1>&);
   2347 
   2348     private:
   2349       template<typename _ForwardIterator,
   2350 	       typename _UniformRandomNumberGenerator>
   2351 	void
   2352 	__generate_impl(_ForwardIterator __f, _ForwardIterator __t,
   2353 			_UniformRandomNumberGenerator& __urng,
   2354 			const param_type& __p);
   2355 
   2356       param_type _M_param;
   2357 
   2358       __gnu_cxx::arcsine_distribution<result_type> _M_ad;
   2359       std::exponential_distribution<result_type> _M_ed;
   2360     };
   2361 
   2362   /**
   2363    * @brief Return true if two Hoyt distributions are not equal.
   2364    */
   2365   template<typename _RealType>
   2366     inline bool
   2367     operator!=(const hoyt_distribution<_RealType>& __d1,
   2368 	       const hoyt_distribution<_RealType>& __d2)
   2369     { return !(__d1 == __d2); }
   2370 
   2371 
   2372   /**
   2373    * @brief A triangular distribution for random numbers.
   2374    *
   2375    * The formula for the triangular probability density function is
   2376    * @f[
   2377    *                  / 0                          for x < a
   2378    *     p(x|a,b,c) = | \frac{2(x-a)}{(c-a)(b-a)}  for a <= x <= b
   2379    *                  | \frac{2(c-x)}{(c-a)(c-b)}  for b < x <= c
   2380    *                  \ 0                          for c < x
   2381    * @f]
   2382    *
   2383    * <table border=1 cellpadding=10 cellspacing=0>
   2384    * <caption align=top>Distribution Statistics</caption>
   2385    * <tr><td>Mean</td><td>@f$ \frac{a+b+c}{2} @f$</td></tr>
   2386    * <tr><td>Variance</td><td>@f$ \frac{a^2+b^2+c^2-ab-ac-bc}
   2387    *                                   {18}@f$</td></tr>
   2388    * <tr><td>Range</td><td>@f$[a, c]@f$</td></tr>
   2389    * </table>
   2390    */
   2391   template<typename _RealType = double>
   2392     class triangular_distribution
   2393     {
   2394       static_assert(std::is_floating_point<_RealType>::value,
   2395 		    "template argument not a floating point type");
   2396 
   2397     public:
   2398       /** The type of the range of the distribution. */
   2399       typedef _RealType result_type;
   2400 
   2401       /** Parameter type. */
   2402       struct param_type
   2403       {
   2404 	friend class triangular_distribution<_RealType>;
   2405 
   2406 	explicit
   2407 	param_type(_RealType __a = _RealType(0),
   2408 		   _RealType __b = _RealType(0.5),
   2409 		   _RealType __c = _RealType(1))
   2410 	: _M_a(__a), _M_b(__b), _M_c(__c)
   2411 	{
   2412 	  __glibcxx_assert(_M_a <= _M_b);
   2413 	  __glibcxx_assert(_M_b <= _M_c);
   2414 	  __glibcxx_assert(_M_a < _M_c);
   2415 
   2416 	  _M_r_ab = (_M_b - _M_a) / (_M_c - _M_a);
   2417 	  _M_f_ab_ac = (_M_b - _M_a) * (_M_c - _M_a);
   2418 	  _M_f_bc_ac = (_M_c - _M_b) * (_M_c - _M_a);
   2419 	}
   2420 
   2421 	_RealType
   2422 	a() const
   2423 	{ return _M_a; }
   2424 
   2425 	_RealType
   2426 	b() const
   2427 	{ return _M_b; }
   2428 
   2429 	_RealType
   2430 	c() const
   2431 	{ return _M_c; }
   2432 
   2433 	friend bool
   2434 	operator==(const param_type& __p1, const param_type& __p2)
   2435 	{
   2436 	  return (__p1._M_a == __p2._M_a && __p1._M_b == __p2._M_b
   2437 		  && __p1._M_c == __p2._M_c);
   2438 	}
   2439 
   2440 	friend bool
   2441 	operator!=(const param_type& __p1, const param_type& __p2)
   2442 	{ return !(__p1 == __p2); }
   2443 
   2444       private:
   2445 
   2446 	_RealType _M_a;
   2447 	_RealType _M_b;
   2448 	_RealType _M_c;
   2449 	_RealType _M_r_ab;
   2450 	_RealType _M_f_ab_ac;
   2451 	_RealType _M_f_bc_ac;
   2452       };
   2453 
   2454       /**
   2455        * @brief Constructs a triangle distribution with parameters
   2456        * @f$ a @f$, @f$ b @f$ and @f$ c @f$.
   2457        */
   2458       explicit
   2459       triangular_distribution(result_type __a = result_type(0),
   2460 			      result_type __b = result_type(0.5),
   2461 			      result_type __c = result_type(1))
   2462       : _M_param(__a, __b, __c)
   2463       { }
   2464 
   2465       explicit
   2466       triangular_distribution(const param_type& __p)
   2467       : _M_param(__p)
   2468       { }
   2469 
   2470       /**
   2471        * @brief Resets the distribution state.
   2472        */
   2473       void
   2474       reset()
   2475       { }
   2476 
   2477       /**
   2478        * @brief Returns the @f$ a @f$ of the distribution.
   2479        */
   2480       result_type
   2481       a() const
   2482       { return _M_param.a(); }
   2483 
   2484       /**
   2485        * @brief Returns the @f$ b @f$ of the distribution.
   2486        */
   2487       result_type
   2488       b() const
   2489       { return _M_param.b(); }
   2490 
   2491       /**
   2492        * @brief Returns the @f$ c @f$ of the distribution.
   2493        */
   2494       result_type
   2495       c() const
   2496       { return _M_param.c(); }
   2497 
   2498       /**
   2499        * @brief Returns the parameter set of the distribution.
   2500        */
   2501       param_type
   2502       param() const
   2503       { return _M_param; }
   2504 
   2505       /**
   2506        * @brief Sets the parameter set of the distribution.
   2507        * @param __param The new parameter set of the distribution.
   2508        */
   2509       void
   2510       param(const param_type& __param)
   2511       { _M_param = __param; }
   2512 
   2513       /**
   2514        * @brief Returns the greatest lower bound value of the distribution.
   2515        */
   2516       result_type
   2517       min() const
   2518       { return _M_param._M_a; }
   2519 
   2520       /**
   2521        * @brief Returns the least upper bound value of the distribution.
   2522        */
   2523       result_type
   2524       max() const
   2525       { return _M_param._M_c; }
   2526 
   2527       /**
   2528        * @brief Generating functions.
   2529        */
   2530       template<typename _UniformRandomNumberGenerator>
   2531 	result_type
   2532 	operator()(_UniformRandomNumberGenerator& __urng)
   2533 	{ return this->operator()(__urng, _M_param); }
   2534 
   2535       template<typename _UniformRandomNumberGenerator>
   2536 	result_type
   2537 	operator()(_UniformRandomNumberGenerator& __urng,
   2538 		   const param_type& __p)
   2539 	{
   2540 	  std::__detail::_Adaptor<_UniformRandomNumberGenerator, result_type>
   2541 	    __aurng(__urng);
   2542 	  result_type __rnd = __aurng();
   2543 	  if (__rnd <= __p._M_r_ab)
   2544 	    return __p.a() + std::sqrt(__rnd * __p._M_f_ab_ac);
   2545 	  else
   2546 	    return __p.c() - std::sqrt((result_type(1) - __rnd)
   2547 				       * __p._M_f_bc_ac);
   2548 	}
   2549 
   2550       template<typename _ForwardIterator,
   2551 	       typename _UniformRandomNumberGenerator>
   2552 	void
   2553 	__generate(_ForwardIterator __f, _ForwardIterator __t,
   2554 		   _UniformRandomNumberGenerator& __urng)
   2555 	{ this->__generate(__f, __t, __urng, _M_param); }
   2556 
   2557       template<typename _ForwardIterator,
   2558 	       typename _UniformRandomNumberGenerator>
   2559 	void
   2560 	__generate(_ForwardIterator __f, _ForwardIterator __t,
   2561 		   _UniformRandomNumberGenerator& __urng,
   2562 		   const param_type& __p)
   2563 	{ this->__generate_impl(__f, __t, __urng, __p); }
   2564 
   2565       template<typename _UniformRandomNumberGenerator>
   2566 	void
   2567 	__generate(result_type* __f, result_type* __t,
   2568 		   _UniformRandomNumberGenerator& __urng,
   2569 		   const param_type& __p)
   2570 	{ this->__generate_impl(__f, __t, __urng, __p); }
   2571 
   2572       /**
   2573        * @brief Return true if two triangle distributions have the same
   2574        *        parameters and the sequences that would be generated
   2575        *        are equal.
   2576        */
   2577       friend bool
   2578       operator==(const triangular_distribution& __d1,
   2579 		 const triangular_distribution& __d2)
   2580       { return __d1._M_param == __d2._M_param; }
   2581 
   2582       /**
   2583        * @brief Inserts a %triangular_distribution random number distribution
   2584        * @p __x into the output stream @p __os.
   2585        *
   2586        * @param __os An output stream.
   2587        * @param __x  A %triangular_distribution random number distribution.
   2588        *
   2589        * @returns The output stream with the state of @p __x inserted or in
   2590        * an error state.
   2591        */
   2592       template<typename _RealType1, typename _CharT, typename _Traits>
   2593 	friend std::basic_ostream<_CharT, _Traits>&
   2594 	operator<<(std::basic_ostream<_CharT, _Traits>& __os,
   2595 		   const __gnu_cxx::triangular_distribution<_RealType1>& __x);
   2596 
   2597       /**
   2598        * @brief Extracts a %triangular_distribution random number distribution
   2599        * @p __x from the input stream @p __is.
   2600        *
   2601        * @param __is An input stream.
   2602        * @param __x  A %triangular_distribution random number generator engine.
   2603        *
   2604        * @returns The input stream with @p __x extracted or in an error state.
   2605        */
   2606       template<typename _RealType1, typename _CharT, typename _Traits>
   2607 	friend std::basic_istream<_CharT, _Traits>&
   2608 	operator>>(std::basic_istream<_CharT, _Traits>& __is,
   2609 		   __gnu_cxx::triangular_distribution<_RealType1>& __x);
   2610 
   2611     private:
   2612       template<typename _ForwardIterator,
   2613 	       typename _UniformRandomNumberGenerator>
   2614 	void
   2615 	__generate_impl(_ForwardIterator __f, _ForwardIterator __t,
   2616 			_UniformRandomNumberGenerator& __urng,
   2617 			const param_type& __p);
   2618 
   2619       param_type _M_param;
   2620     };
   2621 
   2622   /**
   2623    * @brief Return true if two triangle distributions are different.
   2624    */
   2625   template<typename _RealType>
   2626     inline bool
   2627     operator!=(const __gnu_cxx::triangular_distribution<_RealType>& __d1,
   2628 	       const __gnu_cxx::triangular_distribution<_RealType>& __d2)
   2629     { return !(__d1 == __d2); }
   2630 
   2631 
   2632   /**
   2633    * @brief A von Mises distribution for random numbers.
   2634    *
   2635    * The formula for the von Mises probability density function is
   2636    * @f[
   2637    *     p(x|\mu,\kappa) = \frac{e^{\kappa \cos(x-\mu)}}
   2638    *                            {2\pi I_0(\kappa)}
   2639    * @f]
   2640    *
   2641    * The generating functions use the method according to:
   2642    *
   2643    * D. J. Best and N. I. Fisher, 1979. "Efficient Simulation of the
   2644    * von Mises Distribution", Journal of the Royal Statistical Society.
   2645    * Series C (Applied Statistics), Vol. 28, No. 2, pp. 152-157.
   2646    *
   2647    * <table border=1 cellpadding=10 cellspacing=0>
   2648    * <caption align=top>Distribution Statistics</caption>
   2649    * <tr><td>Mean</td><td>@f$ \mu @f$</td></tr>
   2650    * <tr><td>Variance</td><td>@f$ 1-I_1(\kappa)/I_0(\kappa) @f$</td></tr>
   2651    * <tr><td>Range</td><td>@f$[-\pi, \pi]@f$</td></tr>
   2652    * </table>
   2653    */
   2654   template<typename _RealType = double>
   2655     class von_mises_distribution
   2656     {
   2657       static_assert(std::is_floating_point<_RealType>::value,
   2658 		    "template argument not a floating point type");
   2659 
   2660     public:
   2661       /** The type of the range of the distribution. */
   2662       typedef _RealType result_type;
   2663       /** Parameter type. */
   2664       struct param_type
   2665       {
   2666 	friend class von_mises_distribution<_RealType>;
   2667 
   2668 	explicit
   2669 	param_type(_RealType __mu = _RealType(0),
   2670 		   _RealType __kappa = _RealType(1))
   2671 	: _M_mu(__mu), _M_kappa(__kappa)
   2672 	{
   2673 	  const _RealType __pi = __gnu_cxx::__math_constants<_RealType>::__pi;
   2674 	  __glibcxx_assert(_M_mu >= -__pi && _M_mu <= __pi);
   2675 	  __glibcxx_assert(_M_kappa >= _RealType(0));
   2676 
   2677 	  auto __tau = std::sqrt(_RealType(4) * _M_kappa * _M_kappa
   2678 				 + _RealType(1)) + _RealType(1);
   2679 	  auto __rho = ((__tau - std::sqrt(_RealType(2) * __tau))
   2680 			/ (_RealType(2) * _M_kappa));
   2681 	  _M_r = (_RealType(1) + __rho * __rho) / (_RealType(2) * __rho);
   2682 	}
   2683 
   2684 	_RealType
   2685 	mu() const
   2686 	{ return _M_mu; }
   2687 
   2688 	_RealType
   2689 	kappa() const
   2690 	{ return _M_kappa; }
   2691 
   2692 	friend bool
   2693 	operator==(const param_type& __p1, const param_type& __p2)
   2694 	{ return __p1._M_mu == __p2._M_mu && __p1._M_kappa == __p2._M_kappa; }
   2695 
   2696 	friend bool
   2697 	operator!=(const param_type& __p1, const param_type& __p2)
   2698 	{ return !(__p1 == __p2); }
   2699 
   2700       private:
   2701 	_RealType _M_mu;
   2702 	_RealType _M_kappa;
   2703 	_RealType _M_r;
   2704       };
   2705 
   2706       /**
   2707        * @brief Constructs a von Mises distribution with parameters
   2708        * @f$\mu@f$ and @f$\kappa@f$.
   2709        */
   2710       explicit
   2711       von_mises_distribution(result_type __mu = result_type(0),
   2712 			     result_type __kappa = result_type(1))
   2713 	: _M_param(__mu, __kappa)
   2714       { }
   2715 
   2716       explicit
   2717       von_mises_distribution(const param_type& __p)
   2718       : _M_param(__p)
   2719       { }
   2720 
   2721       /**
   2722        * @brief Resets the distribution state.
   2723        */
   2724       void
   2725       reset()
   2726       { }
   2727 
   2728       /**
   2729        * @brief Returns the @f$ \mu @f$ of the distribution.
   2730        */
   2731       result_type
   2732       mu() const
   2733       { return _M_param.mu(); }
   2734 
   2735       /**
   2736        * @brief Returns the @f$ \kappa @f$ of the distribution.
   2737        */
   2738       result_type
   2739       kappa() const
   2740       { return _M_param.kappa(); }
   2741 
   2742       /**
   2743        * @brief Returns the parameter set of the distribution.
   2744        */
   2745       param_type
   2746       param() const
   2747       { return _M_param; }
   2748 
   2749       /**
   2750        * @brief Sets the parameter set of the distribution.
   2751        * @param __param The new parameter set of the distribution.
   2752        */
   2753       void
   2754       param(const param_type& __param)
   2755       { _M_param = __param; }
   2756 
   2757       /**
   2758        * @brief Returns the greatest lower bound value of the distribution.
   2759        */
   2760       result_type
   2761       min() const
   2762       {
   2763 	return -__gnu_cxx::__math_constants<result_type>::__pi;
   2764       }
   2765 
   2766       /**
   2767        * @brief Returns the least upper bound value of the distribution.
   2768        */
   2769       result_type
   2770       max() const
   2771       {
   2772 	return __gnu_cxx::__math_constants<result_type>::__pi;
   2773       }
   2774 
   2775       /**
   2776        * @brief Generating functions.
   2777        */
   2778       template<typename _UniformRandomNumberGenerator>
   2779 	result_type
   2780 	operator()(_UniformRandomNumberGenerator& __urng)
   2781 	{ return this->operator()(__urng, _M_param); }
   2782 
   2783       template<typename _UniformRandomNumberGenerator>
   2784 	result_type
   2785 	operator()(_UniformRandomNumberGenerator& __urng,
   2786 		   const param_type& __p);
   2787 
   2788       template<typename _ForwardIterator,
   2789 	       typename _UniformRandomNumberGenerator>
   2790 	void
   2791 	__generate(_ForwardIterator __f, _ForwardIterator __t,
   2792 		   _UniformRandomNumberGenerator& __urng)
   2793 	{ this->__generate(__f, __t, __urng, _M_param); }
   2794 
   2795       template<typename _ForwardIterator,
   2796 	       typename _UniformRandomNumberGenerator>
   2797 	void
   2798 	__generate(_ForwardIterator __f, _ForwardIterator __t,
   2799 		   _UniformRandomNumberGenerator& __urng,
   2800 		   const param_type& __p)
   2801 	{ this->__generate_impl(__f, __t, __urng, __p); }
   2802 
   2803       template<typename _UniformRandomNumberGenerator>
   2804 	void
   2805 	__generate(result_type* __f, result_type* __t,
   2806 		   _UniformRandomNumberGenerator& __urng,
   2807 		   const param_type& __p)
   2808 	{ this->__generate_impl(__f, __t, __urng, __p); }
   2809 
   2810       /**
   2811        * @brief Return true if two von Mises distributions have the same
   2812        *        parameters and the sequences that would be generated
   2813        *        are equal.
   2814        */
   2815       friend bool
   2816       operator==(const von_mises_distribution& __d1,
   2817 		 const von_mises_distribution& __d2)
   2818       { return __d1._M_param == __d2._M_param; }
   2819 
   2820       /**
   2821        * @brief Inserts a %von_mises_distribution random number distribution
   2822        * @p __x into the output stream @p __os.
   2823        *
   2824        * @param __os An output stream.
   2825        * @param __x  A %von_mises_distribution random number distribution.
   2826        *
   2827        * @returns The output stream with the state of @p __x inserted or in
   2828        * an error state.
   2829        */
   2830       template<typename _RealType1, typename _CharT, typename _Traits>
   2831 	friend std::basic_ostream<_CharT, _Traits>&
   2832 	operator<<(std::basic_ostream<_CharT, _Traits>& __os,
   2833 		   const __gnu_cxx::von_mises_distribution<_RealType1>& __x);
   2834 
   2835       /**
   2836        * @brief Extracts a %von_mises_distribution random number distribution
   2837        * @p __x from the input stream @p __is.
   2838        *
   2839        * @param __is An input stream.
   2840        * @param __x  A %von_mises_distribution random number generator engine.
   2841        *
   2842        * @returns The input stream with @p __x extracted or in an error state.
   2843        */
   2844       template<typename _RealType1, typename _CharT, typename _Traits>
   2845 	friend std::basic_istream<_CharT, _Traits>&
   2846 	operator>>(std::basic_istream<_CharT, _Traits>& __is,
   2847 		   __gnu_cxx::von_mises_distribution<_RealType1>& __x);
   2848 
   2849     private:
   2850       template<typename _ForwardIterator,
   2851 	       typename _UniformRandomNumberGenerator>
   2852 	void
   2853 	__generate_impl(_ForwardIterator __f, _ForwardIterator __t,
   2854 			_UniformRandomNumberGenerator& __urng,
   2855 			const param_type& __p);
   2856 
   2857       param_type _M_param;
   2858     };
   2859 
   2860   /**
   2861    * @brief Return true if two von Mises distributions are different.
   2862    */
   2863   template<typename _RealType>
   2864     inline bool
   2865     operator!=(const __gnu_cxx::von_mises_distribution<_RealType>& __d1,
   2866 	       const __gnu_cxx::von_mises_distribution<_RealType>& __d2)
   2867     { return !(__d1 == __d2); }
   2868 
   2869 
   2870   /**
   2871    * @brief A discrete hypergeometric random number distribution.
   2872    *
   2873    * The hypergeometric distribution is a discrete probability distribution
   2874    * that describes the probability of @p k successes in @p n draws @a without
   2875    * replacement from a finite population of size @p N containing exactly @p K
   2876    * successes.
   2877    *
   2878    * The formula for the hypergeometric probability density function is
   2879    * @f[
   2880    *   p(k|N,K,n) = \frac{\binom{K}{k} \binom{N-K}{n-k}}{\binom{N}{n}}
   2881    * @f]
   2882    * where @f$N@f$ is the total population of the distribution,
   2883    * @f$K@f$ is the total population of the distribution.
   2884    *
   2885    * <table border=1 cellpadding=10 cellspacing=0>
   2886    * <caption align=top>Distribution Statistics</caption>
   2887    * <tr><td>Mean</td><td>@f$ n\frac{K}{N} @f$</td></tr>
   2888    * <tr><td>Variance</td><td>@f$ n\frac{K}{N}\frac{N-K}{N}\frac{N-n}{N-1}
   2889    *   @f$</td></tr>
   2890    * <tr><td>Range</td><td>@f$[max(0, n+K-N), min(K, n)]@f$</td></tr>
   2891    * </table>
   2892    */
   2893   template<typename _UIntType = unsigned int>
   2894     class hypergeometric_distribution
   2895     {
   2896       static_assert(std::is_unsigned<_UIntType>::value, "template argument "
   2897 		    "substituting _UIntType not an unsigned integral type");
   2898 
   2899     public:
   2900       /** The type of the range of the distribution. */
   2901       typedef _UIntType  result_type;
   2902 
   2903       /** Parameter type. */
   2904       struct param_type
   2905       {
   2906 	typedef hypergeometric_distribution<_UIntType> distribution_type;
   2907 	friend class hypergeometric_distribution<_UIntType>;
   2908 
   2909 	explicit
   2910 	param_type(result_type __N = 10, result_type __K = 5,
   2911 		   result_type __n = 1)
   2912 	: _M_N{__N}, _M_K{__K}, _M_n{__n}
   2913 	{
   2914 	  __glibcxx_assert(_M_N >= _M_K);
   2915 	  __glibcxx_assert(_M_N >= _M_n);
   2916 	}
   2917 
   2918 	result_type
   2919 	total_size() const
   2920 	{ return _M_N; }
   2921 
   2922 	result_type
   2923 	successful_size() const
   2924 	{ return _M_K; }
   2925 
   2926 	result_type
   2927 	unsuccessful_size() const
   2928 	{ return _M_N - _M_K; }
   2929 
   2930 	result_type
   2931 	total_draws() const
   2932 	{ return _M_n; }
   2933 
   2934 	friend bool
   2935 	operator==(const param_type& __p1, const param_type& __p2)
   2936 	{ return (__p1._M_N == __p2._M_N)
   2937 	      && (__p1._M_K == __p2._M_K)
   2938 	      && (__p1._M_n == __p2._M_n); }
   2939 
   2940 	friend bool
   2941 	operator!=(const param_type& __p1, const param_type& __p2)
   2942 	{ return !(__p1 == __p2); }
   2943 
   2944       private:
   2945 
   2946 	result_type _M_N;
   2947 	result_type _M_K;
   2948 	result_type _M_n;
   2949       };
   2950 
   2951       // constructors and member function
   2952       explicit
   2953       hypergeometric_distribution(result_type __N = 10, result_type __K = 5,
   2954 				  result_type __n = 1)
   2955       : _M_param{__N, __K, __n}
   2956       { }
   2957 
   2958       explicit
   2959       hypergeometric_distribution(const param_type& __p)
   2960       : _M_param{__p}
   2961       { }
   2962 
   2963       /**
   2964        * @brief Resets the distribution state.
   2965        */
   2966       void
   2967       reset()
   2968       { }
   2969 
   2970       /**
   2971        * @brief Returns the distribution parameter @p N,
   2972        *	the total number of items.
   2973        */
   2974       result_type
   2975       total_size() const
   2976       { return this->_M_param.total_size(); }
   2977 
   2978       /**
   2979        * @brief Returns the distribution parameter @p K,
   2980        *	the total number of successful items.
   2981        */
   2982       result_type
   2983       successful_size() const
   2984       { return this->_M_param.successful_size(); }
   2985 
   2986       /**
   2987        * @brief Returns the total number of unsuccessful items @f$ N - K @f$.
   2988        */
   2989       result_type
   2990       unsuccessful_size() const
   2991       { return this->_M_param.unsuccessful_size(); }
   2992 
   2993       /**
   2994        * @brief Returns the distribution parameter @p n,
   2995        *	the total number of draws.
   2996        */
   2997       result_type
   2998       total_draws() const
   2999       { return this->_M_param.total_draws(); }
   3000 
   3001       /**
   3002        * @brief Returns the parameter set of the distribution.
   3003        */
   3004       param_type
   3005       param() const
   3006       { return this->_M_param; }
   3007 
   3008       /**
   3009        * @brief Sets the parameter set of the distribution.
   3010        * @param __param The new parameter set of the distribution.
   3011        */
   3012       void
   3013       param(const param_type& __param)
   3014       { this->_M_param = __param; }
   3015 
   3016       /**
   3017        * @brief Returns the greatest lower bound value of the distribution.
   3018        */
   3019       result_type
   3020       min() const
   3021       {
   3022 	using _IntType = typename std::make_signed<result_type>::type;
   3023 	return static_cast<result_type>(std::max(static_cast<_IntType>(0),
   3024 				static_cast<_IntType>(this->total_draws()
   3025 						- this->unsuccessful_size())));
   3026       }
   3027 
   3028       /**
   3029        * @brief Returns the least upper bound value of the distribution.
   3030        */
   3031       result_type
   3032       max() const
   3033       { return std::min(this->successful_size(), this->total_draws()); }
   3034 
   3035       /**
   3036        * @brief Generating functions.
   3037        */
   3038       template<typename _UniformRandomNumberGenerator>
   3039 	result_type
   3040 	operator()(_UniformRandomNumberGenerator& __urng)
   3041 	{ return this->operator()(__urng, this->_M_param); }
   3042 
   3043       template<typename _UniformRandomNumberGenerator>
   3044 	result_type
   3045 	operator()(_UniformRandomNumberGenerator& __urng,
   3046 		   const param_type& __p);
   3047 
   3048       template<typename _ForwardIterator,
   3049 	       typename _UniformRandomNumberGenerator>
   3050 	void
   3051 	__generate(_ForwardIterator __f, _ForwardIterator __t,
   3052 		   _UniformRandomNumberGenerator& __urng)
   3053 	{ this->__generate(__f, __t, __urng, this->_M_param); }
   3054 
   3055       template<typename _ForwardIterator,
   3056 	       typename _UniformRandomNumberGenerator>
   3057 	void
   3058 	__generate(_ForwardIterator __f, _ForwardIterator __t,
   3059 		   _UniformRandomNumberGenerator& __urng,
   3060 		   const param_type& __p)
   3061 	{ this->__generate_impl(__f, __t, __urng, __p); }
   3062 
   3063       template<typename _UniformRandomNumberGenerator>
   3064 	void
   3065 	__generate(result_type* __f, result_type* __t,
   3066 		   _UniformRandomNumberGenerator& __urng,
   3067 		   const param_type& __p)
   3068 	{ this->__generate_impl(__f, __t, __urng, __p); }
   3069 
   3070        /**
   3071 	* @brief Return true if two hypergeometric distributions have the same
   3072 	*        parameters and the sequences that would be generated
   3073 	*        are equal.
   3074 	*/
   3075       friend bool
   3076       operator==(const hypergeometric_distribution& __d1,
   3077 		 const hypergeometric_distribution& __d2)
   3078       { return __d1._M_param == __d2._M_param; }
   3079 
   3080       /**
   3081        * @brief Inserts a %hypergeometric_distribution random number
   3082        * distribution @p __x into the output stream @p __os.
   3083        *
   3084        * @param __os An output stream.
   3085        * @param __x  A %hypergeometric_distribution random number
   3086        *             distribution.
   3087        *
   3088        * @returns The output stream with the state of @p __x inserted or in
   3089        * an error state.
   3090        */
   3091       template<typename _UIntType1, typename _CharT, typename _Traits>
   3092 	friend std::basic_ostream<_CharT, _Traits>&
   3093 	operator<<(std::basic_ostream<_CharT, _Traits>& __os,
   3094 		   const __gnu_cxx::hypergeometric_distribution<_UIntType1>&
   3095 		   __x);
   3096 
   3097       /**
   3098        * @brief Extracts a %hypergeometric_distribution random number
   3099        * distribution @p __x from the input stream @p __is.
   3100        *
   3101        * @param __is An input stream.
   3102        * @param __x  A %hypergeometric_distribution random number generator
   3103        *             distribution.
   3104        *
   3105        * @returns The input stream with @p __x extracted or in an error
   3106        *          state.
   3107        */
   3108       template<typename _UIntType1, typename _CharT, typename _Traits>
   3109 	friend std::basic_istream<_CharT, _Traits>&
   3110 	operator>>(std::basic_istream<_CharT, _Traits>& __is,
   3111 		   __gnu_cxx::hypergeometric_distribution<_UIntType1>& __x);
   3112 
   3113     private:
   3114 
   3115       template<typename _ForwardIterator,
   3116 	       typename _UniformRandomNumberGenerator>
   3117 	void
   3118 	__generate_impl(_ForwardIterator __f, _ForwardIterator __t,
   3119 			_UniformRandomNumberGenerator& __urng,
   3120 			const param_type& __p);
   3121 
   3122       param_type _M_param;
   3123     };
   3124 
   3125   /**
   3126    * @brief Return true if two hypergeometric distributions are different.
   3127    */
   3128   template<typename _UIntType>
   3129     inline bool
   3130     operator!=(const __gnu_cxx::hypergeometric_distribution<_UIntType>& __d1,
   3131 	       const __gnu_cxx::hypergeometric_distribution<_UIntType>& __d2)
   3132     { return !(__d1 == __d2); }
   3133 
   3134   /**
   3135    * @brief A logistic continuous distribution for random numbers.
   3136    *
   3137    * The formula for the logistic probability density function is
   3138    * @f[
   3139    *     p(x|\a,\b) = \frac{e^{(x - a)/b}}{b[1 + e^{(x - a)/b}]^2}
   3140    * @f]
   3141    * where @f$b > 0@f$.
   3142    *
   3143    * The formula for the logistic probability function is
   3144    * @f[
   3145    *     cdf(x|\a,\b) = \frac{e^{(x - a)/b}}{1 + e^{(x - a)/b}}
   3146    * @f]
   3147    * where @f$b > 0@f$.
   3148    *
   3149    * <table border=1 cellpadding=10 cellspacing=0>
   3150    * <caption align=top>Distribution Statistics</caption>
   3151    * <tr><td>Mean</td><td>@f$a@f$</td></tr>
   3152    * <tr><td>Variance</td><td>@f$b^2\pi^2/3@f$</td></tr>
   3153    * <tr><td>Range</td><td>@f$[0, \infty)@f$</td></tr>
   3154    * </table>
   3155    */
   3156   template<typename _RealType = double>
   3157     class
   3158     logistic_distribution
   3159     {
   3160       static_assert(std::is_floating_point<_RealType>::value,
   3161 		    "template argument not a floating point type");
   3162 
   3163     public:
   3164       /** The type of the range of the distribution. */
   3165       typedef _RealType result_type;
   3166 
   3167       /** Parameter type. */
   3168       struct param_type
   3169       {
   3170 	typedef logistic_distribution<result_type> distribution_type;
   3171 
   3172 	param_type(result_type __a = result_type(0),
   3173 		   result_type __b = result_type(1))
   3174 	: _M_a(__a), _M_b(__b)
   3175 	{
   3176 	  __glibcxx_assert(_M_b > result_type(0));
   3177 	}
   3178 
   3179 	result_type
   3180 	a() const
   3181 	{ return _M_a; }
   3182 
   3183 	result_type
   3184 	b() const
   3185 	{ return _M_b; }
   3186 
   3187 	friend bool
   3188 	operator==(const param_type& __p1, const param_type& __p2)
   3189 	{ return __p1._M_a == __p2._M_a && __p1._M_b == __p2._M_b; }
   3190 
   3191 	friend bool
   3192 	operator!=(const param_type& __p1, const param_type& __p2)
   3193 	{ return !(__p1 == __p2); }
   3194 
   3195       private:
   3196 	void _M_initialize();
   3197 
   3198 	result_type _M_a;
   3199 	result_type _M_b;
   3200       };
   3201 
   3202       /**
   3203        * @brief Constructors.
   3204        */
   3205       explicit
   3206       logistic_distribution(result_type __a = result_type(0),
   3207 			    result_type __b = result_type(1))
   3208       : _M_param(__a, __b)
   3209       { }
   3210 
   3211       explicit
   3212       logistic_distribution(const param_type& __p)
   3213       : _M_param(__p)
   3214       { }
   3215 
   3216       /**
   3217        * @brief Resets the distribution state.
   3218        */
   3219       void
   3220       reset()
   3221       { }
   3222 
   3223       /**
   3224        * @brief Return the parameters of the distribution.
   3225        */
   3226       result_type
   3227       a() const
   3228       { return _M_param.a(); }
   3229 
   3230       result_type
   3231       b() const
   3232       { return _M_param.b(); }
   3233 
   3234       /**
   3235        * @brief Returns the parameter set of the distribution.
   3236        */
   3237       param_type
   3238       param() const
   3239       { return _M_param; }
   3240 
   3241       /**
   3242        * @brief Sets the parameter set of the distribution.
   3243        * @param __param The new parameter set of the distribution.
   3244        */
   3245       void
   3246       param(const param_type& __param)
   3247       { _M_param = __param; }
   3248 
   3249       /**
   3250        * @brief Returns the greatest lower bound value of the distribution.
   3251        */
   3252       result_type
   3253       min() const
   3254       { return -std::numeric_limits<result_type>::max(); }
   3255 
   3256       /**
   3257        * @brief Returns the least upper bound value of the distribution.
   3258        */
   3259       result_type
   3260       max() const
   3261       { return std::numeric_limits<result_type>::max(); }
   3262 
   3263       /**
   3264        * @brief Generating functions.
   3265        */
   3266       template<typename _UniformRandomNumberGenerator>
   3267 	result_type
   3268 	operator()(_UniformRandomNumberGenerator& __urng)
   3269 	{ return this->operator()(__urng, this->_M_param); }
   3270 
   3271       template<typename _UniformRandomNumberGenerator>
   3272 	result_type
   3273 	operator()(_UniformRandomNumberGenerator&,
   3274 		   const param_type&);
   3275 
   3276       template<typename _ForwardIterator,
   3277 	       typename _UniformRandomNumberGenerator>
   3278 	void
   3279 	__generate(_ForwardIterator __f, _ForwardIterator __t,
   3280 		   _UniformRandomNumberGenerator& __urng)
   3281 	{ this->__generate(__f, __t, __urng, this->param()); }
   3282 
   3283       template<typename _ForwardIterator,
   3284 	       typename _UniformRandomNumberGenerator>
   3285 	void
   3286 	__generate(_ForwardIterator __f, _ForwardIterator __t,
   3287 		   _UniformRandomNumberGenerator& __urng,
   3288 		   const param_type& __p)
   3289 	{ this->__generate_impl(__f, __t, __urng, __p); }
   3290 
   3291       template<typename _UniformRandomNumberGenerator>
   3292 	void
   3293 	__generate(result_type* __f, result_type* __t,
   3294 		   _UniformRandomNumberGenerator& __urng,
   3295 		   const param_type& __p)
   3296 	{ this->__generate_impl(__f, __t, __urng, __p); }
   3297 
   3298       /**
   3299        * @brief Return true if two logistic distributions have
   3300        *        the same parameters and the sequences that would
   3301        *        be generated are equal.
   3302        */
   3303       template<typename _RealType1>
   3304 	friend bool
   3305 	operator==(const logistic_distribution<_RealType1>& __d1,
   3306 		   const logistic_distribution<_RealType1>& __d2)
   3307 	{ return __d1.param() == __d2.param(); }
   3308 
   3309       /**
   3310        * @brief Inserts a %logistic_distribution random number distribution
   3311        * @p __x into the output stream @p __os.
   3312        *
   3313        * @param __os An output stream.
   3314        * @param __x  A %logistic_distribution random number distribution.
   3315        *
   3316        * @returns The output stream with the state of @p __x inserted or in
   3317        * an error state.
   3318        */
   3319       template<typename _RealType1, typename _CharT, typename _Traits>
   3320 	friend std::basic_ostream<_CharT, _Traits>&
   3321 	operator<<(std::basic_ostream<_CharT, _Traits>&,
   3322 		   const logistic_distribution<_RealType1>&);
   3323 
   3324       /**
   3325        * @brief Extracts a %logistic_distribution random number distribution
   3326        * @p __x from the input stream @p __is.
   3327        *
   3328        * @param __is An input stream.
   3329        * @param __x A %logistic_distribution random number
   3330        *            generator engine.
   3331        *
   3332        * @returns The input stream with @p __x extracted or in an error state.
   3333        */
   3334       template<typename _RealType1, typename _CharT, typename _Traits>
   3335 	friend std::basic_istream<_CharT, _Traits>&
   3336 	operator>>(std::basic_istream<_CharT, _Traits>&,
   3337 		   logistic_distribution<_RealType1>&);
   3338 
   3339     private:
   3340       template<typename _ForwardIterator,
   3341 	       typename _UniformRandomNumberGenerator>
   3342 	void
   3343 	__generate_impl(_ForwardIterator __f, _ForwardIterator __t,
   3344 			_UniformRandomNumberGenerator& __urng,
   3345 			const param_type& __p);
   3346 
   3347       param_type _M_param;
   3348     };
   3349 
   3350   /**
   3351    * @brief Return true if two logistic distributions are not equal.
   3352    */
   3353   template<typename _RealType1>
   3354     inline bool
   3355     operator!=(const logistic_distribution<_RealType1>& __d1,
   3356 	       const logistic_distribution<_RealType1>& __d2)
   3357     { return !(__d1 == __d2); }
   3358 
   3359 
   3360   /**
   3361    * @brief A distribution for random coordinates on a unit sphere.
   3362    *
   3363    * The method used in the generation function is attributed by Donald Knuth
   3364    * to G. W. Brown, Modern Mathematics for the Engineer (1956).
   3365    */
   3366   template<std::size_t _Dimen, typename _RealType = double>
   3367     class uniform_on_sphere_distribution
   3368     {
   3369       static_assert(std::is_floating_point<_RealType>::value,
   3370 		    "template argument not a floating point type");
   3371       static_assert(_Dimen != 0, "dimension is zero");
   3372 
   3373     public:
   3374       /** The type of the range of the distribution. */
   3375       typedef std::array<_RealType, _Dimen> result_type;
   3376 
   3377       /** Parameter type. */
   3378       struct param_type
   3379       {
   3380 	explicit
   3381 	param_type()
   3382 	{ }
   3383 
   3384 	friend bool
   3385 	operator==(const param_type&, const param_type&)
   3386 	{ return true; }
   3387 
   3388 	friend bool
   3389 	operator!=(const param_type&, const param_type&)
   3390 	{ return false; }
   3391       };
   3392 
   3393       /**
   3394        * @brief Constructs a uniform on sphere distribution.
   3395        */
   3396       explicit
   3397       uniform_on_sphere_distribution()
   3398       : _M_param(), _M_nd()
   3399       { }
   3400 
   3401       explicit
   3402       uniform_on_sphere_distribution(const param_type& __p)
   3403       : _M_param(__p), _M_nd()
   3404       { }
   3405 
   3406       /**
   3407        * @brief Resets the distribution state.
   3408        */
   3409       void
   3410       reset()
   3411       { _M_nd.reset(); }
   3412 
   3413       /**
   3414        * @brief Returns the parameter set of the distribution.
   3415        */
   3416       param_type
   3417       param() const
   3418       { return _M_param; }
   3419 
   3420       /**
   3421        * @brief Sets the parameter set of the distribution.
   3422        * @param __param The new parameter set of the distribution.
   3423        */
   3424       void
   3425       param(const param_type& __param)
   3426       { _M_param = __param; }
   3427 
   3428       /**
   3429        * @brief Returns the greatest lower bound value of the distribution.
   3430        * This function makes no sense for this distribution.
   3431        */
   3432       result_type
   3433       min() const
   3434       {
   3435 	result_type __res;
   3436 	__res.fill(0);
   3437 	return __res;
   3438       }
   3439 
   3440       /**
   3441        * @brief Returns the least upper bound value of the distribution.
   3442        * This function makes no sense for this distribution.
   3443        */
   3444       result_type
   3445       max() const
   3446       {
   3447 	result_type __res;
   3448 	__res.fill(0);
   3449 	return __res;
   3450       }
   3451 
   3452       /**
   3453        * @brief Generating functions.
   3454        */
   3455       template<typename _UniformRandomNumberGenerator>
   3456 	result_type
   3457 	operator()(_UniformRandomNumberGenerator& __urng)
   3458 	{ return this->operator()(__urng, _M_param); }
   3459 
   3460       template<typename _UniformRandomNumberGenerator>
   3461 	result_type
   3462 	operator()(_UniformRandomNumberGenerator& __urng,
   3463 		   const param_type& __p);
   3464 
   3465       template<typename _ForwardIterator,
   3466 	       typename _UniformRandomNumberGenerator>
   3467 	void
   3468 	__generate(_ForwardIterator __f, _ForwardIterator __t,
   3469 		   _UniformRandomNumberGenerator& __urng)
   3470 	{ this->__generate(__f, __t, __urng, this->param()); }
   3471 
   3472       template<typename _ForwardIterator,
   3473 	       typename _UniformRandomNumberGenerator>
   3474 	void
   3475 	__generate(_ForwardIterator __f, _ForwardIterator __t,
   3476 		   _UniformRandomNumberGenerator& __urng,
   3477 		   const param_type& __p)
   3478 	{ this->__generate_impl(__f, __t, __urng, __p); }
   3479 
   3480       template<typename _UniformRandomNumberGenerator>
   3481 	void
   3482 	__generate(result_type* __f, result_type* __t,
   3483 		   _UniformRandomNumberGenerator& __urng,
   3484 		   const param_type& __p)
   3485 	{ this->__generate_impl(__f, __t, __urng, __p); }
   3486 
   3487       /**
   3488        * @brief Return true if two uniform on sphere distributions have
   3489        *        the same parameters and the sequences that would be
   3490        *        generated are equal.
   3491        */
   3492       friend bool
   3493       operator==(const uniform_on_sphere_distribution& __d1,
   3494 		 const uniform_on_sphere_distribution& __d2)
   3495       { return __d1._M_nd == __d2._M_nd; }
   3496 
   3497       /**
   3498        * @brief Inserts a %uniform_on_sphere_distribution random number
   3499        *        distribution @p __x into the output stream @p __os.
   3500        *
   3501        * @param __os An output stream.
   3502        * @param __x  A %uniform_on_sphere_distribution random number
   3503        *             distribution.
   3504        *
   3505        * @returns The output stream with the state of @p __x inserted or in
   3506        * an error state.
   3507        */
   3508       template<size_t _Dimen1, typename _RealType1, typename _CharT,
   3509 	       typename _Traits>
   3510 	friend std::basic_ostream<_CharT, _Traits>&
   3511 	operator<<(std::basic_ostream<_CharT, _Traits>& __os,
   3512 		   const __gnu_cxx::uniform_on_sphere_distribution<_Dimen1,
   3513 								   _RealType1>&
   3514 		   __x);
   3515 
   3516       /**
   3517        * @brief Extracts a %uniform_on_sphere_distribution random number
   3518        *        distribution
   3519        * @p __x from the input stream @p __is.
   3520        *
   3521        * @param __is An input stream.
   3522        * @param __x  A %uniform_on_sphere_distribution random number
   3523        *             generator engine.
   3524        *
   3525        * @returns The input stream with @p __x extracted or in an error state.
   3526        */
   3527       template<std::size_t _Dimen1, typename _RealType1, typename _CharT,
   3528 	       typename _Traits>
   3529 	friend std::basic_istream<_CharT, _Traits>&
   3530 	operator>>(std::basic_istream<_CharT, _Traits>& __is,
   3531 		   __gnu_cxx::uniform_on_sphere_distribution<_Dimen1,
   3532 							     _RealType1>& __x);
   3533 
   3534     private:
   3535       template<typename _ForwardIterator,
   3536 	       typename _UniformRandomNumberGenerator>
   3537 	void
   3538 	__generate_impl(_ForwardIterator __f, _ForwardIterator __t,
   3539 			_UniformRandomNumberGenerator& __urng,
   3540 			const param_type& __p);
   3541 
   3542       param_type _M_param;
   3543       std::normal_distribution<_RealType> _M_nd;
   3544     };
   3545 
   3546   /**
   3547    * @brief Return true if two uniform on sphere distributions are different.
   3548    */
   3549   template<std::size_t _Dimen, typename _RealType>
   3550     inline bool
   3551     operator!=(const __gnu_cxx::uniform_on_sphere_distribution<_Dimen,
   3552 	       _RealType>& __d1,
   3553 	       const __gnu_cxx::uniform_on_sphere_distribution<_Dimen,
   3554 	       _RealType>& __d2)
   3555     { return !(__d1 == __d2); }
   3556 
   3557 
   3558   /**
   3559    * @brief A distribution for random coordinates inside a unit sphere.
   3560    */
   3561   template<std::size_t _Dimen, typename _RealType = double>
   3562     class uniform_inside_sphere_distribution
   3563     {
   3564       static_assert(std::is_floating_point<_RealType>::value,
   3565 		    "template argument not a floating point type");
   3566       static_assert(_Dimen != 0, "dimension is zero");
   3567 
   3568     public:
   3569       /** The type of the range of the distribution. */
   3570       using result_type = std::array<_RealType, _Dimen>;
   3571 
   3572       /** Parameter type. */
   3573       struct param_type
   3574       {
   3575 	using distribution_type
   3576 	  = uniform_inside_sphere_distribution<_Dimen, _RealType>;
   3577 	friend class uniform_inside_sphere_distribution<_Dimen, _RealType>;
   3578 
   3579 	explicit
   3580 	param_type(_RealType __radius = _RealType(1))
   3581 	: _M_radius(__radius)
   3582 	{
   3583 	  __glibcxx_assert(_M_radius > _RealType(0));
   3584 	}
   3585 
   3586 	_RealType
   3587 	radius() const
   3588 	{ return _M_radius; }
   3589 
   3590 	friend bool
   3591 	operator==(const param_type& __p1, const param_type& __p2)
   3592 	{ return __p1._M_radius == __p2._M_radius; }
   3593 
   3594 	friend bool
   3595 	operator!=(const param_type& __p1, const param_type& __p2)
   3596 	{ return !(__p1 == __p2); }
   3597 
   3598       private:
   3599 	_RealType _M_radius;
   3600       };
   3601 
   3602       /**
   3603        * @brief Constructors.
   3604        */
   3605       explicit
   3606       uniform_inside_sphere_distribution(_RealType __radius = _RealType(1))
   3607       : _M_param(__radius), _M_uosd()
   3608       { }
   3609 
   3610       explicit
   3611       uniform_inside_sphere_distribution(const param_type& __p)
   3612       : _M_param(__p), _M_uosd()
   3613       { }
   3614 
   3615       /**
   3616        * @brief Resets the distribution state.
   3617        */
   3618       void
   3619       reset()
   3620       { _M_uosd.reset(); }
   3621 
   3622       /**
   3623        * @brief Returns the @f$radius@f$ of the distribution.
   3624        */
   3625       _RealType
   3626       radius() const
   3627       { return _M_param.radius(); }
   3628 
   3629       /**
   3630        * @brief Returns the parameter set of the distribution.
   3631        */
   3632       param_type
   3633       param() const
   3634       { return _M_param; }
   3635 
   3636       /**
   3637        * @brief Sets the parameter set of the distribution.
   3638        * @param __param The new parameter set of the distribution.
   3639        */
   3640       void
   3641       param(const param_type& __param)
   3642       { _M_param = __param; }
   3643 
   3644       /**
   3645        * @brief Returns the greatest lower bound value of the distribution.
   3646        * This function makes no sense for this distribution.
   3647        */
   3648       result_type
   3649       min() const
   3650       {
   3651 	result_type __res;
   3652 	__res.fill(0);
   3653 	return __res;
   3654       }
   3655 
   3656       /**
   3657        * @brief Returns the least upper bound value of the distribution.
   3658        * This function makes no sense for this distribution.
   3659        */
   3660       result_type
   3661       max() const
   3662       {
   3663 	result_type __res;
   3664 	__res.fill(0);
   3665 	return __res;
   3666       }
   3667 
   3668       /**
   3669        * @brief Generating functions.
   3670        */
   3671       template<typename _UniformRandomNumberGenerator>
   3672 	result_type
   3673 	operator()(_UniformRandomNumberGenerator& __urng)
   3674 	{ return this->operator()(__urng, _M_param); }
   3675 
   3676       template<typename _UniformRandomNumberGenerator>
   3677 	result_type
   3678 	operator()(_UniformRandomNumberGenerator& __urng,
   3679 		   const param_type& __p);
   3680 
   3681       template<typename _ForwardIterator,
   3682 	       typename _UniformRandomNumberGenerator>
   3683 	void
   3684 	__generate(_ForwardIterator __f, _ForwardIterator __t,
   3685 		   _UniformRandomNumberGenerator& __urng)
   3686 	{ this->__generate(__f, __t, __urng, this->param()); }
   3687 
   3688       template<typename _ForwardIterator,
   3689 	       typename _UniformRandomNumberGenerator>
   3690 	void
   3691 	__generate(_ForwardIterator __f, _ForwardIterator __t,
   3692 		   _UniformRandomNumberGenerator& __urng,
   3693 		   const param_type& __p)
   3694 	{ this->__generate_impl(__f, __t, __urng, __p); }
   3695 
   3696       template<typename _UniformRandomNumberGenerator>
   3697 	void
   3698 	__generate(result_type* __f, result_type* __t,
   3699 		   _UniformRandomNumberGenerator& __urng,
   3700 		   const param_type& __p)
   3701 	{ this->__generate_impl(__f, __t, __urng, __p); }
   3702 
   3703       /**
   3704        * @brief Return true if two uniform on sphere distributions have
   3705        *        the same parameters and the sequences that would be
   3706        *        generated are equal.
   3707        */
   3708       friend bool
   3709       operator==(const uniform_inside_sphere_distribution& __d1,
   3710 		 const uniform_inside_sphere_distribution& __d2)
   3711       { return __d1._M_param == __d2._M_param && __d1._M_uosd == __d2._M_uosd; }
   3712 
   3713       /**
   3714        * @brief Inserts a %uniform_inside_sphere_distribution random number
   3715        *        distribution @p __x into the output stream @p __os.
   3716        *
   3717        * @param __os An output stream.
   3718        * @param __x  A %uniform_inside_sphere_distribution random number
   3719        *             distribution.
   3720        *
   3721        * @returns The output stream with the state of @p __x inserted or in
   3722        * an error state.
   3723        */
   3724       template<size_t _Dimen1, typename _RealType1, typename _CharT,
   3725 	       typename _Traits>
   3726 	friend std::basic_ostream<_CharT, _Traits>&
   3727 	operator<<(std::basic_ostream<_CharT, _Traits>& __os,
   3728 		   const __gnu_cxx::uniform_inside_sphere_distribution<_Dimen1,
   3729 								   _RealType1>&
   3730 		   );
   3731 
   3732       /**
   3733        * @brief Extracts a %uniform_inside_sphere_distribution random number
   3734        *        distribution
   3735        * @p __x from the input stream @p __is.
   3736        *
   3737        * @param __is An input stream.
   3738        * @param __x  A %uniform_inside_sphere_distribution random number
   3739        *             generator engine.
   3740        *
   3741        * @returns The input stream with @p __x extracted or in an error state.
   3742        */
   3743       template<std::size_t _Dimen1, typename _RealType1, typename _CharT,
   3744 	       typename _Traits>
   3745 	friend std::basic_istream<_CharT, _Traits>&
   3746 	operator>>(std::basic_istream<_CharT, _Traits>& __is,
   3747 		   __gnu_cxx::uniform_inside_sphere_distribution<_Dimen1,
   3748 								 _RealType1>&);
   3749 
   3750     private:
   3751       template<typename _ForwardIterator,
   3752 	       typename _UniformRandomNumberGenerator>
   3753 	void
   3754 	__generate_impl(_ForwardIterator __f, _ForwardIterator __t,
   3755 			_UniformRandomNumberGenerator& __urng,
   3756 			const param_type& __p);
   3757 
   3758       param_type _M_param;
   3759       uniform_on_sphere_distribution<_Dimen, _RealType> _M_uosd;
   3760     };
   3761 
   3762   /**
   3763    * @brief Return true if two uniform on sphere distributions are different.
   3764    */
   3765   template<std::size_t _Dimen, typename _RealType>
   3766     inline bool
   3767     operator!=(const __gnu_cxx::uniform_inside_sphere_distribution<_Dimen,
   3768 	       _RealType>& __d1,
   3769 	       const __gnu_cxx::uniform_inside_sphere_distribution<_Dimen,
   3770 	       _RealType>& __d2)
   3771     { return !(__d1 == __d2); }
   3772 
   3773 _GLIBCXX_END_NAMESPACE_VERSION
   3774 } // namespace __gnu_cxx
   3775 
   3776 #include "ext/opt_random.h"
   3777 #include "random.tcc"
   3778 
   3779 #endif // _GLIBCXX_USE_C99_STDINT_TR1 && UINT32_C
   3780 
   3781 #endif // C++11
   3782 
   3783 #endif // _EXT_RANDOM
   3784