1// Random number extensions -*- C++ -*- 2 3// Copyright (C) 2012-2016 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 48namespace __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 /** Parameter type. */ 408 struct param_type 409 { 410 typedef beta_distribution<_RealType> distribution_type; 411 friend class beta_distribution<_RealType>; 412 413 explicit 414 param_type(_RealType __alpha_val = _RealType(1), 415 _RealType __beta_val = _RealType(1)) 416 : _M_alpha(__alpha_val), _M_beta(__beta_val) 417 { 418 __glibcxx_assert(_M_alpha > _RealType(0)); 419 __glibcxx_assert(_M_beta > _RealType(0)); 420 } 421 422 _RealType 423 alpha() const 424 { return _M_alpha; } 425 426 _RealType 427 beta() const 428 { return _M_beta; } 429 430 friend bool 431 operator==(const param_type& __p1, const param_type& __p2) 432 { return (__p1._M_alpha == __p2._M_alpha 433 && __p1._M_beta == __p2._M_beta); } 434 435 private: 436 void 437 _M_initialize(); 438 439 _RealType _M_alpha; 440 _RealType _M_beta; 441 }; 442 443 public: 444 /** 445 * @brief Constructs a beta distribution with parameters 446 * @f$\alpha@f$ and @f$\beta@f$. 447 */ 448 explicit 449 beta_distribution(_RealType __alpha_val = _RealType(1), 450 _RealType __beta_val = _RealType(1)) 451 : _M_param(__alpha_val, __beta_val) 452 { } 453 454 explicit 455 beta_distribution(const param_type& __p) 456 : _M_param(__p) 457 { } 458 459 /** 460 * @brief Resets the distribution state. 461 */ 462 void 463 reset() 464 { } 465 466 /** 467 * @brief Returns the @f$\alpha@f$ of the distribution. 468 */ 469 _RealType 470 alpha() const 471 { return _M_param.alpha(); } 472 473 /** 474 * @brief Returns the @f$\beta@f$ of the distribution. 475 */ 476 _RealType 477 beta() const 478 { return _M_param.beta(); } 479 480 /** 481 * @brief Returns the parameter set of the distribution. 482 */ 483 param_type 484 param() const 485 { return _M_param; } 486 487 /** 488 * @brief Sets the parameter set of the distribution. 489 * @param __param The new parameter set of the distribution. 490 */ 491 void 492 param(const param_type& __param) 493 { _M_param = __param; } 494 495 /** 496 * @brief Returns the greatest lower bound value of the distribution. 497 */ 498 result_type 499 min() const 500 { return result_type(0); } 501 502 /** 503 * @brief Returns the least upper bound value of the distribution. 504 */ 505 result_type 506 max() const 507 { return result_type(1); } 508 509 /** 510 * @brief Generating functions. 511 */ 512 template<typename _UniformRandomNumberGenerator> 513 result_type 514 operator()(_UniformRandomNumberGenerator& __urng) 515 { return this->operator()(__urng, _M_param); } 516 517 template<typename _UniformRandomNumberGenerator> 518 result_type 519 operator()(_UniformRandomNumberGenerator& __urng, 520 const param_type& __p); 521 522 template<typename _ForwardIterator, 523 typename _UniformRandomNumberGenerator> 524 void 525 __generate(_ForwardIterator __f, _ForwardIterator __t, 526 _UniformRandomNumberGenerator& __urng) 527 { this->__generate(__f, __t, __urng, _M_param); } 528 529 template<typename _ForwardIterator, 530 typename _UniformRandomNumberGenerator> 531 void 532 __generate(_ForwardIterator __f, _ForwardIterator __t, 533 _UniformRandomNumberGenerator& __urng, 534 const param_type& __p) 535 { this->__generate_impl(__f, __t, __urng, __p); } 536 537 template<typename _UniformRandomNumberGenerator> 538 void 539 __generate(result_type* __f, result_type* __t, 540 _UniformRandomNumberGenerator& __urng, 541 const param_type& __p) 542 { this->__generate_impl(__f, __t, __urng, __p); } 543 544 /** 545 * @brief Return true if two beta distributions have the same 546 * parameters and the sequences that would be generated 547 * are equal. 548 */ 549 friend bool 550 operator==(const beta_distribution& __d1, 551 const beta_distribution& __d2) 552 { return __d1._M_param == __d2._M_param; } 553 554 /** 555 * @brief Inserts a %beta_distribution random number distribution 556 * @p __x into the output stream @p __os. 557 * 558 * @param __os An output stream. 559 * @param __x A %beta_distribution random number distribution. 560 * 561 * @returns The output stream with the state of @p __x inserted or in 562 * an error state. 563 */ 564 template<typename _RealType1, typename _CharT, typename _Traits> 565 friend std::basic_ostream<_CharT, _Traits>& 566 operator<<(std::basic_ostream<_CharT, _Traits>& __os, 567 const __gnu_cxx::beta_distribution<_RealType1>& __x); 568 569 /** 570 * @brief Extracts a %beta_distribution random number distribution 571 * @p __x from the input stream @p __is. 572 * 573 * @param __is An input stream. 574 * @param __x A %beta_distribution random number generator engine. 575 * 576 * @returns The input stream with @p __x extracted or in an error state. 577 */ 578 template<typename _RealType1, typename _CharT, typename _Traits> 579 friend std::basic_istream<_CharT, _Traits>& 580 operator>>(std::basic_istream<_CharT, _Traits>& __is, 581 __gnu_cxx::beta_distribution<_RealType1>& __x); 582 583 private: 584 template<typename _ForwardIterator, 585 typename _UniformRandomNumberGenerator> 586 void 587 __generate_impl(_ForwardIterator __f, _ForwardIterator __t, 588 _UniformRandomNumberGenerator& __urng, 589 const param_type& __p); 590 591 param_type _M_param; 592 }; 593 594 /** 595 * @brief Return true if two beta distributions are different. 596 */ 597 template<typename _RealType> 598 inline bool 599 operator!=(const __gnu_cxx::beta_distribution<_RealType>& __d1, 600 const __gnu_cxx::beta_distribution<_RealType>& __d2) 601 { return !(__d1 == __d2); } 602 603 604 /** 605 * @brief A multi-variate normal continuous distribution for random numbers. 606 * 607 * The formula for the normal probability density function is 608 * @f[ 609 * p(\overrightarrow{x}|\overrightarrow{\mu },\Sigma) = 610 * \frac{1}{\sqrt{(2\pi )^k\det(\Sigma))}} 611 * e^{-\frac{1}{2}(\overrightarrow{x}-\overrightarrow{\mu})^\text{T} 612 * \Sigma ^{-1}(\overrightarrow{x}-\overrightarrow{\mu})} 613 * @f] 614 * 615 * where @f$\overrightarrow{x}@f$ and @f$\overrightarrow{\mu}@f$ are 616 * vectors of dimension @f$k@f$ and @f$\Sigma@f$ is the covariance 617 * matrix (which must be positive-definite). 618 */ 619 template<std::size_t _Dimen, typename _RealType = double> 620 class normal_mv_distribution 621 { 622 static_assert(std::is_floating_point<_RealType>::value, 623 "template argument not a floating point type"); 624 static_assert(_Dimen != 0, "dimension is zero"); 625 626 public: 627 /** The type of the range of the distribution. */ 628 typedef std::array<_RealType, _Dimen> result_type; 629 /** Parameter type. */ 630 class param_type 631 { 632 static constexpr size_t _M_t_size = _Dimen * (_Dimen + 1) / 2; 633 634 public: 635 typedef normal_mv_distribution<_Dimen, _RealType> distribution_type; 636 friend class normal_mv_distribution<_Dimen, _RealType>; 637 638 param_type() 639 { 640 std::fill(_M_mean.begin(), _M_mean.end(), _RealType(0)); 641 auto __it = _M_t.begin(); 642 for (size_t __i = 0; __i < _Dimen; ++__i) 643 { 644 std::fill_n(__it, __i, _RealType(0)); 645 __it += __i; 646 *__it++ = _RealType(1); 647 } 648 } 649 650 template<typename _ForwardIterator1, typename _ForwardIterator2> 651 param_type(_ForwardIterator1 __meanbegin, 652 _ForwardIterator1 __meanend, 653 _ForwardIterator2 __varcovbegin, 654 _ForwardIterator2 __varcovend) 655 { 656 __glibcxx_function_requires(_ForwardIteratorConcept< 657 _ForwardIterator1>) 658 __glibcxx_function_requires(_ForwardIteratorConcept< 659 _ForwardIterator2>) 660 _GLIBCXX_DEBUG_ASSERT(std::distance(__meanbegin, __meanend) 661 <= _Dimen); 662 const auto __dist = std::distance(__varcovbegin, __varcovend); 663 _GLIBCXX_DEBUG_ASSERT(__dist == _Dimen * _Dimen 664 || __dist == _Dimen * (_Dimen + 1) / 2 665 || __dist == _Dimen); 666 667 if (__dist == _Dimen * _Dimen) 668 _M_init_full(__meanbegin, __meanend, __varcovbegin, __varcovend); 669 else if (__dist == _Dimen * (_Dimen + 1) / 2) 670 _M_init_lower(__meanbegin, __meanend, __varcovbegin, __varcovend); 671 else 672 { 673 __glibcxx_assert(__dist == _Dimen); 674 _M_init_diagonal(__meanbegin, __meanend, 675 __varcovbegin, __varcovend); 676 } 677 } 678 679 param_type(std::initializer_list<_RealType> __mean, 680 std::initializer_list<_RealType> __varcov) 681 { 682 _GLIBCXX_DEBUG_ASSERT(__mean.size() <= _Dimen); 683 _GLIBCXX_DEBUG_ASSERT(__varcov.size() == _Dimen * _Dimen 684 || __varcov.size() == _Dimen * (_Dimen + 1) / 2 685 || __varcov.size() == _Dimen); 686 687 if (__varcov.size() == _Dimen * _Dimen) 688 _M_init_full(__mean.begin(), __mean.end(), 689 __varcov.begin(), __varcov.end()); 690 else if (__varcov.size() == _Dimen * (_Dimen + 1) / 2) 691 _M_init_lower(__mean.begin(), __mean.end(), 692 __varcov.begin(), __varcov.end()); 693 else 694 { 695 __glibcxx_assert(__varcov.size() == _Dimen); 696 _M_init_diagonal(__mean.begin(), __mean.end(), 697 __varcov.begin(), __varcov.end()); 698 } 699 } 700 701 std::array<_RealType, _Dimen> 702 mean() const 703 { return _M_mean; } 704 705 std::array<_RealType, _M_t_size> 706 varcov() const 707 { return _M_t; } 708 709 friend bool 710 operator==(const param_type& __p1, const param_type& __p2) 711 { return __p1._M_mean == __p2._M_mean && __p1._M_t == __p2._M_t; } 712 713 private: 714 template <typename _InputIterator1, typename _InputIterator2> 715 void _M_init_full(_InputIterator1 __meanbegin, 716 _InputIterator1 __meanend, 717 _InputIterator2 __varcovbegin, 718 _InputIterator2 __varcovend); 719 template <typename _InputIterator1, typename _InputIterator2> 720 void _M_init_lower(_InputIterator1 __meanbegin, 721 _InputIterator1 __meanend, 722 _InputIterator2 __varcovbegin, 723 _InputIterator2 __varcovend); 724 template <typename _InputIterator1, typename _InputIterator2> 725 void _M_init_diagonal(_InputIterator1 __meanbegin, 726 _InputIterator1 __meanend, 727 _InputIterator2 __varbegin, 728 _InputIterator2 __varend); 729 730 std::array<_RealType, _Dimen> _M_mean; 731 std::array<_RealType, _M_t_size> _M_t; 732 }; 733 734 public: 735 normal_mv_distribution() 736 : _M_param(), _M_nd() 737 { } 738 739 template<typename _ForwardIterator1, typename _ForwardIterator2> 740 normal_mv_distribution(_ForwardIterator1 __meanbegin, 741 _ForwardIterator1 __meanend, 742 _ForwardIterator2 __varcovbegin, 743 _ForwardIterator2 __varcovend) 744 : _M_param(__meanbegin, __meanend, __varcovbegin, __varcovend), 745 _M_nd() 746 { } 747 748 normal_mv_distribution(std::initializer_list<_RealType> __mean, 749 std::initializer_list<_RealType> __varcov) 750 : _M_param(__mean, __varcov), _M_nd() 751 { } 752 753 explicit 754 normal_mv_distribution(const param_type& __p) 755 : _M_param(__p), _M_nd() 756 { } 757 758 /** 759 * @brief Resets the distribution state. 760 */ 761 void 762 reset() 763 { _M_nd.reset(); } 764 765 /** 766 * @brief Returns the mean of the distribution. 767 */ 768 result_type 769 mean() const 770 { return _M_param.mean(); } 771 772 /** 773 * @brief Returns the compact form of the variance/covariance 774 * matrix of the distribution. 775 */ 776 std::array<_RealType, _Dimen * (_Dimen + 1) / 2> 777 varcov() const 778 { return _M_param.varcov(); } 779 780 /** 781 * @brief Returns the parameter set of the distribution. 782 */ 783 param_type 784 param() const 785 { return _M_param; } 786 787 /** 788 * @brief Sets the parameter set of the distribution. 789 * @param __param The new parameter set of the distribution. 790 */ 791 void 792 param(const param_type& __param) 793 { _M_param = __param; } 794 795 /** 796 * @brief Returns the greatest lower bound value of the distribution. 797 */ 798 result_type 799 min() const 800 { result_type __res; 801 __res.fill(std::numeric_limits<_RealType>::lowest()); 802 return __res; } 803 804 /** 805 * @brief Returns the least upper bound value of the distribution. 806 */ 807 result_type 808 max() const 809 { result_type __res; 810 __res.fill(std::numeric_limits<_RealType>::max()); 811 return __res; } 812 813 /** 814 * @brief Generating functions. 815 */ 816 template<typename _UniformRandomNumberGenerator> 817 result_type 818 operator()(_UniformRandomNumberGenerator& __urng) 819 { return this->operator()(__urng, _M_param); } 820 821 template<typename _UniformRandomNumberGenerator> 822 result_type 823 operator()(_UniformRandomNumberGenerator& __urng, 824 const param_type& __p); 825 826 template<typename _ForwardIterator, 827 typename _UniformRandomNumberGenerator> 828 void 829 __generate(_ForwardIterator __f, _ForwardIterator __t, 830 _UniformRandomNumberGenerator& __urng) 831 { return this->__generate_impl(__f, __t, __urng, _M_param); } 832 833 template<typename _ForwardIterator, 834 typename _UniformRandomNumberGenerator> 835 void 836 __generate(_ForwardIterator __f, _ForwardIterator __t, 837 _UniformRandomNumberGenerator& __urng, 838 const param_type& __p) 839 { return this->__generate_impl(__f, __t, __urng, __p); } 840 841 /** 842 * @brief Return true if two multi-variant normal distributions have 843 * the same parameters and the sequences that would 844 * be generated are equal. 845 */ 846 template<size_t _Dimen1, typename _RealType1> 847 friend bool 848 operator==(const 849 __gnu_cxx::normal_mv_distribution<_Dimen1, _RealType1>& 850 __d1, 851 const 852 __gnu_cxx::normal_mv_distribution<_Dimen1, _RealType1>& 853 __d2); 854 855 /** 856 * @brief Inserts a %normal_mv_distribution random number distribution 857 * @p __x into the output stream @p __os. 858 * 859 * @param __os An output stream. 860 * @param __x A %normal_mv_distribution random number distribution. 861 * 862 * @returns The output stream with the state of @p __x inserted or in 863 * an error state. 864 */ 865 template<size_t _Dimen1, typename _RealType1, 866 typename _CharT, typename _Traits> 867 friend std::basic_ostream<_CharT, _Traits>& 868 operator<<(std::basic_ostream<_CharT, _Traits>& __os, 869 const 870 __gnu_cxx::normal_mv_distribution<_Dimen1, _RealType1>& 871 __x); 872 873 /** 874 * @brief Extracts a %normal_mv_distribution random number distribution 875 * @p __x from the input stream @p __is. 876 * 877 * @param __is An input stream. 878 * @param __x A %normal_mv_distribution random number generator engine. 879 * 880 * @returns The input stream with @p __x extracted or in an error 881 * state. 882 */ 883 template<size_t _Dimen1, typename _RealType1, 884 typename _CharT, typename _Traits> 885 friend std::basic_istream<_CharT, _Traits>& 886 operator>>(std::basic_istream<_CharT, _Traits>& __is, 887 __gnu_cxx::normal_mv_distribution<_Dimen1, _RealType1>& 888 __x); 889 890 private: 891 template<typename _ForwardIterator, 892 typename _UniformRandomNumberGenerator> 893 void 894 __generate_impl(_ForwardIterator __f, _ForwardIterator __t, 895 _UniformRandomNumberGenerator& __urng, 896 const param_type& __p); 897 898 param_type _M_param; 899 std::normal_distribution<_RealType> _M_nd; 900 }; 901 902 /** 903 * @brief Return true if two multi-variate normal distributions are 904 * different. 905 */ 906 template<size_t _Dimen, typename _RealType> 907 inline bool 908 operator!=(const __gnu_cxx::normal_mv_distribution<_Dimen, _RealType>& 909 __d1, 910 const __gnu_cxx::normal_mv_distribution<_Dimen, _RealType>& 911 __d2) 912 { return !(__d1 == __d2); } 913 914 915 /** 916 * @brief A Rice continuous distribution for random numbers. 917 * 918 * The formula for the Rice probability density function is 919 * @f[ 920 * p(x|\nu,\sigma) = \frac{x}{\sigma^2} 921 * \exp\left(-\frac{x^2+\nu^2}{2\sigma^2}\right) 922 * I_0\left(\frac{x \nu}{\sigma^2}\right) 923 * @f] 924 * where @f$I_0(z)@f$ is the modified Bessel function of the first kind 925 * of order 0 and @f$\nu >= 0@f$ and @f$\sigma > 0@f$. 926 * 927 * <table border=1 cellpadding=10 cellspacing=0> 928 * <caption align=top>Distribution Statistics</caption> 929 * <tr><td>Mean</td><td>@f$\sqrt{\pi/2}L_{1/2}(-\nu^2/2\sigma^2)@f$</td></tr> 930 * <tr><td>Variance</td><td>@f$2\sigma^2 + \nu^2 931 * + (\pi\sigma^2/2)L^2_{1/2}(-\nu^2/2\sigma^2)@f$</td></tr> 932 * <tr><td>Range</td><td>@f$[0, \infty)@f$</td></tr> 933 * </table> 934 * where @f$L_{1/2}(x)@f$ is the Laguerre polynomial of order 1/2. 935 */ 936 template<typename _RealType = double> 937 class 938 rice_distribution 939 { 940 static_assert(std::is_floating_point<_RealType>::value, 941 "template argument not a floating point type"); 942 public: 943 /** The type of the range of the distribution. */ 944 typedef _RealType result_type; 945 /** Parameter type. */ 946 struct param_type 947 { 948 typedef rice_distribution<result_type> distribution_type; 949 950 param_type(result_type __nu_val = result_type(0), 951 result_type __sigma_val = result_type(1)) 952 : _M_nu(__nu_val), _M_sigma(__sigma_val) 953 { 954 __glibcxx_assert(_M_nu >= result_type(0)); 955 __glibcxx_assert(_M_sigma > result_type(0)); 956 } 957 958 result_type 959 nu() const 960 { return _M_nu; } 961 962 result_type 963 sigma() const 964 { return _M_sigma; } 965 966 friend bool 967 operator==(const param_type& __p1, const param_type& __p2) 968 { return __p1._M_nu == __p2._M_nu 969 && __p1._M_sigma == __p2._M_sigma; } 970 971 private: 972 void _M_initialize(); 973 974 result_type _M_nu; 975 result_type _M_sigma; 976 }; 977 978 /** 979 * @brief Constructors. 980 */ 981 explicit 982 rice_distribution(result_type __nu_val = result_type(0), 983 result_type __sigma_val = result_type(1)) 984 : _M_param(__nu_val, __sigma_val), 985 _M_ndx(__nu_val, __sigma_val), 986 _M_ndy(result_type(0), __sigma_val) 987 { } 988 989 explicit 990 rice_distribution(const param_type& __p) 991 : _M_param(__p), 992 _M_ndx(__p.nu(), __p.sigma()), 993 _M_ndy(result_type(0), __p.sigma()) 994 { } 995 996 /** 997 * @brief Resets the distribution state. 998 */ 999 void 1000 reset() 1001 { 1002 _M_ndx.reset(); 1003 _M_ndy.reset(); 1004 } 1005 1006 /** 1007 * @brief Return the parameters of the distribution. 1008 */ 1009 result_type 1010 nu() const 1011 { return _M_param.nu(); } 1012 1013 result_type 1014 sigma() const 1015 { return _M_param.sigma(); } 1016 1017 /** 1018 * @brief Returns the parameter set of the distribution. 1019 */ 1020 param_type 1021 param() const 1022 { return _M_param; } 1023 1024 /** 1025 * @brief Sets the parameter set of the distribution. 1026 * @param __param The new parameter set of the distribution. 1027 */ 1028 void 1029 param(const param_type& __param) 1030 { _M_param = __param; } 1031 1032 /** 1033 * @brief Returns the greatest lower bound value of the distribution. 1034 */ 1035 result_type 1036 min() const 1037 { return result_type(0); } 1038 1039 /** 1040 * @brief Returns the least upper bound value of the distribution. 1041 */ 1042 result_type 1043 max() const 1044 { return std::numeric_limits<result_type>::max(); } 1045 1046 /** 1047 * @brief Generating functions. 1048 */ 1049 template<typename _UniformRandomNumberGenerator> 1050 result_type 1051 operator()(_UniformRandomNumberGenerator& __urng) 1052 { 1053 result_type __x = this->_M_ndx(__urng); 1054 result_type __y = this->_M_ndy(__urng); 1055#if _GLIBCXX_USE_C99_MATH_TR1 1056 return std::hypot(__x, __y); 1057#else 1058 return std::sqrt(__x * __x + __y * __y); 1059#endif 1060 } 1061 1062 template<typename _UniformRandomNumberGenerator> 1063 result_type 1064 operator()(_UniformRandomNumberGenerator& __urng, 1065 const param_type& __p) 1066 { 1067 typename std::normal_distribution<result_type>::param_type 1068 __px(__p.nu(), __p.sigma()), __py(result_type(0), __p.sigma()); 1069 result_type __x = this->_M_ndx(__px, __urng); 1070 result_type __y = this->_M_ndy(__py, __urng); 1071#if _GLIBCXX_USE_C99_MATH_TR1 1072 return std::hypot(__x, __y); 1073#else 1074 return std::sqrt(__x * __x + __y * __y); 1075#endif 1076 } 1077 1078 template<typename _ForwardIterator, 1079 typename _UniformRandomNumberGenerator> 1080 void 1081 __generate(_ForwardIterator __f, _ForwardIterator __t, 1082 _UniformRandomNumberGenerator& __urng) 1083 { this->__generate(__f, __t, __urng, _M_param); } 1084 1085 template<typename _ForwardIterator, 1086 typename _UniformRandomNumberGenerator> 1087 void 1088 __generate(_ForwardIterator __f, _ForwardIterator __t, 1089 _UniformRandomNumberGenerator& __urng, 1090 const param_type& __p) 1091 { this->__generate_impl(__f, __t, __urng, __p); } 1092 1093 template<typename _UniformRandomNumberGenerator> 1094 void 1095 __generate(result_type* __f, result_type* __t, 1096 _UniformRandomNumberGenerator& __urng, 1097 const param_type& __p) 1098 { this->__generate_impl(__f, __t, __urng, __p); } 1099 1100 /** 1101 * @brief Return true if two Rice distributions have 1102 * the same parameters and the sequences that would 1103 * be generated are equal. 1104 */ 1105 friend bool 1106 operator==(const rice_distribution& __d1, 1107 const rice_distribution& __d2) 1108 { return (__d1._M_param == __d2._M_param 1109 && __d1._M_ndx == __d2._M_ndx 1110 && __d1._M_ndy == __d2._M_ndy); } 1111 1112 /** 1113 * @brief Inserts a %rice_distribution random number distribution 1114 * @p __x into the output stream @p __os. 1115 * 1116 * @param __os An output stream. 1117 * @param __x A %rice_distribution random number distribution. 1118 * 1119 * @returns The output stream with the state of @p __x inserted or in 1120 * an error state. 1121 */ 1122 template<typename _RealType1, typename _CharT, typename _Traits> 1123 friend std::basic_ostream<_CharT, _Traits>& 1124 operator<<(std::basic_ostream<_CharT, _Traits>&, 1125 const rice_distribution<_RealType1>&); 1126 1127 /** 1128 * @brief Extracts a %rice_distribution random number distribution 1129 * @p __x from the input stream @p __is. 1130 * 1131 * @param __is An input stream. 1132 * @param __x A %rice_distribution random number 1133 * generator engine. 1134 * 1135 * @returns The input stream with @p __x extracted or in an error state. 1136 */ 1137 template<typename _RealType1, typename _CharT, typename _Traits> 1138 friend std::basic_istream<_CharT, _Traits>& 1139 operator>>(std::basic_istream<_CharT, _Traits>&, 1140 rice_distribution<_RealType1>&); 1141 1142 private: 1143 template<typename _ForwardIterator, 1144 typename _UniformRandomNumberGenerator> 1145 void 1146 __generate_impl(_ForwardIterator __f, _ForwardIterator __t, 1147 _UniformRandomNumberGenerator& __urng, 1148 const param_type& __p); 1149 1150 param_type _M_param; 1151 1152 std::normal_distribution<result_type> _M_ndx; 1153 std::normal_distribution<result_type> _M_ndy; 1154 }; 1155 1156 /** 1157 * @brief Return true if two Rice distributions are not equal. 1158 */ 1159 template<typename _RealType1> 1160 inline bool 1161 operator!=(const rice_distribution<_RealType1>& __d1, 1162 const rice_distribution<_RealType1>& __d2) 1163 { return !(__d1 == __d2); } 1164 1165 1166 /** 1167 * @brief A Nakagami continuous distribution for random numbers. 1168 * 1169 * The formula for the Nakagami probability density function is 1170 * @f[ 1171 * p(x|\mu,\omega) = \frac{2\mu^\mu}{\Gamma(\mu)\omega^\mu} 1172 * x^{2\mu-1}e^{-\mu x / \omega} 1173 * @f] 1174 * where @f$\Gamma(z)@f$ is the gamma function and @f$\mu >= 0.5@f$ 1175 * and @f$\omega > 0@f$. 1176 */ 1177 template<typename _RealType = double> 1178 class 1179 nakagami_distribution 1180 { 1181 static_assert(std::is_floating_point<_RealType>::value, 1182 "template argument not a floating point type"); 1183 1184 public: 1185 /** The type of the range of the distribution. */ 1186 typedef _RealType result_type; 1187 /** Parameter type. */ 1188 struct param_type 1189 { 1190 typedef nakagami_distribution<result_type> distribution_type; 1191 1192 param_type(result_type __mu_val = result_type(1), 1193 result_type __omega_val = result_type(1)) 1194 : _M_mu(__mu_val), _M_omega(__omega_val) 1195 { 1196 __glibcxx_assert(_M_mu >= result_type(0.5L)); 1197 __glibcxx_assert(_M_omega > result_type(0)); 1198 } 1199 1200 result_type 1201 mu() const 1202 { return _M_mu; } 1203 1204 result_type 1205 omega() const 1206 { return _M_omega; } 1207 1208 friend bool 1209 operator==(const param_type& __p1, const param_type& __p2) 1210 { return __p1._M_mu == __p2._M_mu 1211 && __p1._M_omega == __p2._M_omega; } 1212 1213 private: 1214 void _M_initialize(); 1215 1216 result_type _M_mu; 1217 result_type _M_omega; 1218 }; 1219 1220 /** 1221 * @brief Constructors. 1222 */ 1223 explicit 1224 nakagami_distribution(result_type __mu_val = result_type(1), 1225 result_type __omega_val = result_type(1)) 1226 : _M_param(__mu_val, __omega_val), 1227 _M_gd(__mu_val, __omega_val / __mu_val) 1228 { } 1229 1230 explicit 1231 nakagami_distribution(const param_type& __p) 1232 : _M_param(__p), 1233 _M_gd(__p.mu(), __p.omega() / __p.mu()) 1234 { } 1235 1236 /** 1237 * @brief Resets the distribution state. 1238 */ 1239 void 1240 reset() 1241 { _M_gd.reset(); } 1242 1243 /** 1244 * @brief Return the parameters of the distribution. 1245 */ 1246 result_type 1247 mu() const 1248 { return _M_param.mu(); } 1249 1250 result_type 1251 omega() const 1252 { return _M_param.omega(); } 1253 1254 /** 1255 * @brief Returns the parameter set of the distribution. 1256 */ 1257 param_type 1258 param() const 1259 { return _M_param; } 1260 1261 /** 1262 * @brief Sets the parameter set of the distribution. 1263 * @param __param The new parameter set of the distribution. 1264 */ 1265 void 1266 param(const param_type& __param) 1267 { _M_param = __param; } 1268 1269 /** 1270 * @brief Returns the greatest lower bound value of the distribution. 1271 */ 1272 result_type 1273 min() const 1274 { return result_type(0); } 1275 1276 /** 1277 * @brief Returns the least upper bound value of the distribution. 1278 */ 1279 result_type 1280 max() const 1281 { return std::numeric_limits<result_type>::max(); } 1282 1283 /** 1284 * @brief Generating functions. 1285 */ 1286 template<typename _UniformRandomNumberGenerator> 1287 result_type 1288 operator()(_UniformRandomNumberGenerator& __urng) 1289 { return std::sqrt(this->_M_gd(__urng)); } 1290 1291 template<typename _UniformRandomNumberGenerator> 1292 result_type 1293 operator()(_UniformRandomNumberGenerator& __urng, 1294 const param_type& __p) 1295 { 1296 typename std::gamma_distribution<result_type>::param_type 1297 __pg(__p.mu(), __p.omega() / __p.mu()); 1298 return std::sqrt(this->_M_gd(__pg, __urng)); 1299 } 1300 1301 template<typename _ForwardIterator, 1302 typename _UniformRandomNumberGenerator> 1303 void 1304 __generate(_ForwardIterator __f, _ForwardIterator __t, 1305 _UniformRandomNumberGenerator& __urng) 1306 { this->__generate(__f, __t, __urng, _M_param); } 1307 1308 template<typename _ForwardIterator, 1309 typename _UniformRandomNumberGenerator> 1310 void 1311 __generate(_ForwardIterator __f, _ForwardIterator __t, 1312 _UniformRandomNumberGenerator& __urng, 1313 const param_type& __p) 1314 { this->__generate_impl(__f, __t, __urng, __p); } 1315 1316 template<typename _UniformRandomNumberGenerator> 1317 void 1318 __generate(result_type* __f, result_type* __t, 1319 _UniformRandomNumberGenerator& __urng, 1320 const param_type& __p) 1321 { this->__generate_impl(__f, __t, __urng, __p); } 1322 1323 /** 1324 * @brief Return true if two Nakagami distributions have 1325 * the same parameters and the sequences that would 1326 * be generated are equal. 1327 */ 1328 friend bool 1329 operator==(const nakagami_distribution& __d1, 1330 const nakagami_distribution& __d2) 1331 { return (__d1._M_param == __d2._M_param 1332 && __d1._M_gd == __d2._M_gd); } 1333 1334 /** 1335 * @brief Inserts a %nakagami_distribution random number distribution 1336 * @p __x into the output stream @p __os. 1337 * 1338 * @param __os An output stream. 1339 * @param __x A %nakagami_distribution random number distribution. 1340 * 1341 * @returns The output stream with the state of @p __x inserted or in 1342 * an error state. 1343 */ 1344 template<typename _RealType1, typename _CharT, typename _Traits> 1345 friend std::basic_ostream<_CharT, _Traits>& 1346 operator<<(std::basic_ostream<_CharT, _Traits>&, 1347 const nakagami_distribution<_RealType1>&); 1348 1349 /** 1350 * @brief Extracts a %nakagami_distribution random number distribution 1351 * @p __x from the input stream @p __is. 1352 * 1353 * @param __is An input stream. 1354 * @param __x A %nakagami_distribution random number 1355 * generator engine. 1356 * 1357 * @returns The input stream with @p __x extracted or in an error state. 1358 */ 1359 template<typename _RealType1, typename _CharT, typename _Traits> 1360 friend std::basic_istream<_CharT, _Traits>& 1361 operator>>(std::basic_istream<_CharT, _Traits>&, 1362 nakagami_distribution<_RealType1>&); 1363 1364 private: 1365 template<typename _ForwardIterator, 1366 typename _UniformRandomNumberGenerator> 1367 void 1368 __generate_impl(_ForwardIterator __f, _ForwardIterator __t, 1369 _UniformRandomNumberGenerator& __urng, 1370 const param_type& __p); 1371 1372 param_type _M_param; 1373 1374 std::gamma_distribution<result_type> _M_gd; 1375 }; 1376 1377 /** 1378 * @brief Return true if two Nakagami distributions are not equal. 1379 */ 1380 template<typename _RealType> 1381 inline bool 1382 operator!=(const nakagami_distribution<_RealType>& __d1, 1383 const nakagami_distribution<_RealType>& __d2) 1384 { return !(__d1 == __d2); } 1385 1386 1387 /** 1388 * @brief A Pareto continuous distribution for random numbers. 1389 * 1390 * The formula for the Pareto cumulative probability function is 1391 * @f[ 1392 * P(x|\alpha,\mu) = 1 - \left(\frac{\mu}{x}\right)^\alpha 1393 * @f] 1394 * The formula for the Pareto probability density function is 1395 * @f[ 1396 * p(x|\alpha,\mu) = \frac{\alpha + 1}{\mu} 1397 * \left(\frac{\mu}{x}\right)^{\alpha + 1} 1398 * @f] 1399 * where @f$x >= \mu@f$ and @f$\mu > 0@f$, @f$\alpha > 0@f$. 1400 * 1401 * <table border=1 cellpadding=10 cellspacing=0> 1402 * <caption align=top>Distribution Statistics</caption> 1403 * <tr><td>Mean</td><td>@f$\alpha \mu / (\alpha - 1)@f$ 1404 * for @f$\alpha > 1@f$</td></tr> 1405 * <tr><td>Variance</td><td>@f$\alpha \mu^2 / [(\alpha - 1)^2(\alpha - 2)]@f$ 1406 * for @f$\alpha > 2@f$</td></tr> 1407 * <tr><td>Range</td><td>@f$[\mu, \infty)@f$</td></tr> 1408 * </table> 1409 */ 1410 template<typename _RealType = double> 1411 class 1412 pareto_distribution 1413 { 1414 static_assert(std::is_floating_point<_RealType>::value, 1415 "template argument not a floating point type"); 1416 1417 public: 1418 /** The type of the range of the distribution. */ 1419 typedef _RealType result_type; 1420 /** Parameter type. */ 1421 struct param_type 1422 { 1423 typedef pareto_distribution<result_type> distribution_type; 1424 1425 param_type(result_type __alpha_val = result_type(1), 1426 result_type __mu_val = result_type(1)) 1427 : _M_alpha(__alpha_val), _M_mu(__mu_val) 1428 { 1429 __glibcxx_assert(_M_alpha > result_type(0)); 1430 __glibcxx_assert(_M_mu > result_type(0)); 1431 } 1432 1433 result_type 1434 alpha() const 1435 { return _M_alpha; } 1436 1437 result_type 1438 mu() const 1439 { return _M_mu; } 1440 1441 friend bool 1442 operator==(const param_type& __p1, const param_type& __p2) 1443 { return __p1._M_alpha == __p2._M_alpha && __p1._M_mu == __p2._M_mu; } 1444 1445 private: 1446 void _M_initialize(); 1447 1448 result_type _M_alpha; 1449 result_type _M_mu; 1450 }; 1451 1452 /** 1453 * @brief Constructors. 1454 */ 1455 explicit 1456 pareto_distribution(result_type __alpha_val = result_type(1), 1457 result_type __mu_val = result_type(1)) 1458 : _M_param(__alpha_val, __mu_val), 1459 _M_ud() 1460 { } 1461 1462 explicit 1463 pareto_distribution(const param_type& __p) 1464 : _M_param(__p), 1465 _M_ud() 1466 { } 1467 1468 /** 1469 * @brief Resets the distribution state. 1470 */ 1471 void 1472 reset() 1473 { 1474 _M_ud.reset(); 1475 } 1476 1477 /** 1478 * @brief Return the parameters of the distribution. 1479 */ 1480 result_type 1481 alpha() const 1482 { return _M_param.alpha(); } 1483 1484 result_type 1485 mu() const 1486 { return _M_param.mu(); } 1487 1488 /** 1489 * @brief Returns the parameter set of the distribution. 1490 */ 1491 param_type 1492 param() const 1493 { return _M_param; } 1494 1495 /** 1496 * @brief Sets the parameter set of the distribution. 1497 * @param __param The new parameter set of the distribution. 1498 */ 1499 void 1500 param(const param_type& __param) 1501 { _M_param = __param; } 1502 1503 /** 1504 * @brief Returns the greatest lower bound value of the distribution. 1505 */ 1506 result_type 1507 min() const 1508 { return this->mu(); } 1509 1510 /** 1511 * @brief Returns the least upper bound value of the distribution. 1512 */ 1513 result_type 1514 max() const 1515 { return std::numeric_limits<result_type>::max(); } 1516 1517 /** 1518 * @brief Generating functions. 1519 */ 1520 template<typename _UniformRandomNumberGenerator> 1521 result_type 1522 operator()(_UniformRandomNumberGenerator& __urng) 1523 { 1524 return this->mu() * std::pow(this->_M_ud(__urng), 1525 -result_type(1) / this->alpha()); 1526 } 1527 1528 template<typename _UniformRandomNumberGenerator> 1529 result_type 1530 operator()(_UniformRandomNumberGenerator& __urng, 1531 const param_type& __p) 1532 { 1533 return __p.mu() * std::pow(this->_M_ud(__urng), 1534 -result_type(1) / __p.alpha()); 1535 } 1536 1537 template<typename _ForwardIterator, 1538 typename _UniformRandomNumberGenerator> 1539 void 1540 __generate(_ForwardIterator __f, _ForwardIterator __t, 1541 _UniformRandomNumberGenerator& __urng) 1542 { this->__generate(__f, __t, __urng, _M_param); } 1543 1544 template<typename _ForwardIterator, 1545 typename _UniformRandomNumberGenerator> 1546 void 1547 __generate(_ForwardIterator __f, _ForwardIterator __t, 1548 _UniformRandomNumberGenerator& __urng, 1549 const param_type& __p) 1550 { this->__generate_impl(__f, __t, __urng, __p); } 1551 1552 template<typename _UniformRandomNumberGenerator> 1553 void 1554 __generate(result_type* __f, result_type* __t, 1555 _UniformRandomNumberGenerator& __urng, 1556 const param_type& __p) 1557 { this->__generate_impl(__f, __t, __urng, __p); } 1558 1559 /** 1560 * @brief Return true if two Pareto distributions have 1561 * the same parameters and the sequences that would 1562 * be generated are equal. 1563 */ 1564 friend bool 1565 operator==(const pareto_distribution& __d1, 1566 const pareto_distribution& __d2) 1567 { return (__d1._M_param == __d2._M_param 1568 && __d1._M_ud == __d2._M_ud); } 1569 1570 /** 1571 * @brief Inserts a %pareto_distribution random number distribution 1572 * @p __x into the output stream @p __os. 1573 * 1574 * @param __os An output stream. 1575 * @param __x A %pareto_distribution random number distribution. 1576 * 1577 * @returns The output stream with the state of @p __x inserted or in 1578 * an error state. 1579 */ 1580 template<typename _RealType1, typename _CharT, typename _Traits> 1581 friend std::basic_ostream<_CharT, _Traits>& 1582 operator<<(std::basic_ostream<_CharT, _Traits>&, 1583 const pareto_distribution<_RealType1>&); 1584 1585 /** 1586 * @brief Extracts a %pareto_distribution random number distribution 1587 * @p __x from the input stream @p __is. 1588 * 1589 * @param __is An input stream. 1590 * @param __x A %pareto_distribution random number 1591 * generator engine. 1592 * 1593 * @returns The input stream with @p __x extracted or in an error state. 1594 */ 1595 template<typename _RealType1, typename _CharT, typename _Traits> 1596 friend std::basic_istream<_CharT, _Traits>& 1597 operator>>(std::basic_istream<_CharT, _Traits>&, 1598 pareto_distribution<_RealType1>&); 1599 1600 private: 1601 template<typename _ForwardIterator, 1602 typename _UniformRandomNumberGenerator> 1603 void 1604 __generate_impl(_ForwardIterator __f, _ForwardIterator __t, 1605 _UniformRandomNumberGenerator& __urng, 1606 const param_type& __p); 1607 1608 param_type _M_param; 1609 1610 std::uniform_real_distribution<result_type> _M_ud; 1611 }; 1612 1613 /** 1614 * @brief Return true if two Pareto distributions are not equal. 1615 */ 1616 template<typename _RealType> 1617 inline bool 1618 operator!=(const pareto_distribution<_RealType>& __d1, 1619 const pareto_distribution<_RealType>& __d2) 1620 { return !(__d1 == __d2); } 1621 1622 1623 /** 1624 * @brief A K continuous distribution for random numbers. 1625 * 1626 * The formula for the K probability density function is 1627 * @f[ 1628 * p(x|\lambda, \mu, \nu) = \frac{2}{x} 1629 * \left(\frac{\lambda\nu x}{\mu}\right)^{\frac{\lambda + \nu}{2}} 1630 * \frac{1}{\Gamma(\lambda)\Gamma(\nu)} 1631 * K_{\nu - \lambda}\left(2\sqrt{\frac{\lambda\nu x}{\mu}}\right) 1632 * @f] 1633 * where @f$I_0(z)@f$ is the modified Bessel function of the second kind 1634 * of order @f$\nu - \lambda@f$ and @f$\lambda > 0@f$, @f$\mu > 0@f$ 1635 * and @f$\nu > 0@f$. 1636 * 1637 * <table border=1 cellpadding=10 cellspacing=0> 1638 * <caption align=top>Distribution Statistics</caption> 1639 * <tr><td>Mean</td><td>@f$\mu@f$</td></tr> 1640 * <tr><td>Variance</td><td>@f$\mu^2\frac{\lambda + \nu + 1}{\lambda\nu}@f$</td></tr> 1641 * <tr><td>Range</td><td>@f$[0, \infty)@f$</td></tr> 1642 * </table> 1643 */ 1644 template<typename _RealType = double> 1645 class 1646 k_distribution 1647 { 1648 static_assert(std::is_floating_point<_RealType>::value, 1649 "template argument not a floating point type"); 1650 1651 public: 1652 /** The type of the range of the distribution. */ 1653 typedef _RealType result_type; 1654 /** Parameter type. */ 1655 struct param_type 1656 { 1657 typedef k_distribution<result_type> distribution_type; 1658 1659 param_type(result_type __lambda_val = result_type(1), 1660 result_type __mu_val = result_type(1), 1661 result_type __nu_val = result_type(1)) 1662 : _M_lambda(__lambda_val), _M_mu(__mu_val), _M_nu(__nu_val) 1663 { 1664 __glibcxx_assert(_M_lambda > result_type(0)); 1665 __glibcxx_assert(_M_mu > result_type(0)); 1666 __glibcxx_assert(_M_nu > result_type(0)); 1667 } 1668 1669 result_type 1670 lambda() const 1671 { return _M_lambda; } 1672 1673 result_type 1674 mu() const 1675 { return _M_mu; } 1676 1677 result_type 1678 nu() const 1679 { return _M_nu; } 1680 1681 friend bool 1682 operator==(const param_type& __p1, const param_type& __p2) 1683 { return __p1._M_lambda == __p2._M_lambda 1684 && __p1._M_mu == __p2._M_mu 1685 && __p1._M_nu == __p2._M_nu; } 1686 1687 private: 1688 void _M_initialize(); 1689 1690 result_type _M_lambda; 1691 result_type _M_mu; 1692 result_type _M_nu; 1693 }; 1694 1695 /** 1696 * @brief Constructors. 1697 */ 1698 explicit 1699 k_distribution(result_type __lambda_val = result_type(1), 1700 result_type __mu_val = result_type(1), 1701 result_type __nu_val = result_type(1)) 1702 : _M_param(__lambda_val, __mu_val, __nu_val), 1703 _M_gd1(__lambda_val, result_type(1) / __lambda_val), 1704 _M_gd2(__nu_val, __mu_val / __nu_val) 1705 { } 1706 1707 explicit 1708 k_distribution(const param_type& __p) 1709 : _M_param(__p), 1710 _M_gd1(__p.lambda(), result_type(1) / __p.lambda()), 1711 _M_gd2(__p.nu(), __p.mu() / __p.nu()) 1712 { } 1713 1714 /** 1715 * @brief Resets the distribution state. 1716 */ 1717 void 1718 reset() 1719 { 1720 _M_gd1.reset(); 1721 _M_gd2.reset(); 1722 } 1723 1724 /** 1725 * @brief Return the parameters of the distribution. 1726 */ 1727 result_type 1728 lambda() const 1729 { return _M_param.lambda(); } 1730 1731 result_type 1732 mu() const 1733 { return _M_param.mu(); } 1734 1735 result_type 1736 nu() const 1737 { return _M_param.nu(); } 1738 1739 /** 1740 * @brief Returns the parameter set of the distribution. 1741 */ 1742 param_type 1743 param() const 1744 { return _M_param; } 1745 1746 /** 1747 * @brief Sets the parameter set of the distribution. 1748 * @param __param The new parameter set of the distribution. 1749 */ 1750 void 1751 param(const param_type& __param) 1752 { _M_param = __param; } 1753 1754 /** 1755 * @brief Returns the greatest lower bound value of the distribution. 1756 */ 1757 result_type 1758 min() const 1759 { return result_type(0); } 1760 1761 /** 1762 * @brief Returns the least upper bound value of the distribution. 1763 */ 1764 result_type 1765 max() const 1766 { return std::numeric_limits<result_type>::max(); } 1767 1768 /** 1769 * @brief Generating functions. 1770 */ 1771 template<typename _UniformRandomNumberGenerator> 1772 result_type 1773 operator()(_UniformRandomNumberGenerator&); 1774 1775 template<typename _UniformRandomNumberGenerator> 1776 result_type 1777 operator()(_UniformRandomNumberGenerator&, const param_type&); 1778 1779 template<typename _ForwardIterator, 1780 typename _UniformRandomNumberGenerator> 1781 void 1782 __generate(_ForwardIterator __f, _ForwardIterator __t, 1783 _UniformRandomNumberGenerator& __urng) 1784 { this->__generate(__f, __t, __urng, _M_param); } 1785 1786 template<typename _ForwardIterator, 1787 typename _UniformRandomNumberGenerator> 1788 void 1789 __generate(_ForwardIterator __f, _ForwardIterator __t, 1790 _UniformRandomNumberGenerator& __urng, 1791 const param_type& __p) 1792 { this->__generate_impl(__f, __t, __urng, __p); } 1793 1794 template<typename _UniformRandomNumberGenerator> 1795 void 1796 __generate(result_type* __f, result_type* __t, 1797 _UniformRandomNumberGenerator& __urng, 1798 const param_type& __p) 1799 { this->__generate_impl(__f, __t, __urng, __p); } 1800 1801 /** 1802 * @brief Return true if two K distributions have 1803 * the same parameters and the sequences that would 1804 * be generated are equal. 1805 */ 1806 friend bool 1807 operator==(const k_distribution& __d1, 1808 const k_distribution& __d2) 1809 { return (__d1._M_param == __d2._M_param 1810 && __d1._M_gd1 == __d2._M_gd1 1811 && __d1._M_gd2 == __d2._M_gd2); } 1812 1813 /** 1814 * @brief Inserts a %k_distribution random number distribution 1815 * @p __x into the output stream @p __os. 1816 * 1817 * @param __os An output stream. 1818 * @param __x A %k_distribution random number distribution. 1819 * 1820 * @returns The output stream with the state of @p __x inserted or in 1821 * an error state. 1822 */ 1823 template<typename _RealType1, typename _CharT, typename _Traits> 1824 friend std::basic_ostream<_CharT, _Traits>& 1825 operator<<(std::basic_ostream<_CharT, _Traits>&, 1826 const k_distribution<_RealType1>&); 1827 1828 /** 1829 * @brief Extracts a %k_distribution random number distribution 1830 * @p __x from the input stream @p __is. 1831 * 1832 * @param __is An input stream. 1833 * @param __x A %k_distribution random number 1834 * generator engine. 1835 * 1836 * @returns The input stream with @p __x extracted or in an error state. 1837 */ 1838 template<typename _RealType1, typename _CharT, typename _Traits> 1839 friend std::basic_istream<_CharT, _Traits>& 1840 operator>>(std::basic_istream<_CharT, _Traits>&, 1841 k_distribution<_RealType1>&); 1842 1843 private: 1844 template<typename _ForwardIterator, 1845 typename _UniformRandomNumberGenerator> 1846 void 1847 __generate_impl(_ForwardIterator __f, _ForwardIterator __t, 1848 _UniformRandomNumberGenerator& __urng, 1849 const param_type& __p); 1850 1851 param_type _M_param; 1852 1853 std::gamma_distribution<result_type> _M_gd1; 1854 std::gamma_distribution<result_type> _M_gd2; 1855 }; 1856 1857 /** 1858 * @brief Return true if two K distributions are not equal. 1859 */ 1860 template<typename _RealType> 1861 inline bool 1862 operator!=(const k_distribution<_RealType>& __d1, 1863 const k_distribution<_RealType>& __d2) 1864 { return !(__d1 == __d2); } 1865 1866 1867 /** 1868 * @brief An arcsine continuous distribution for random numbers. 1869 * 1870 * The formula for the arcsine probability density function is 1871 * @f[ 1872 * p(x|a,b) = \frac{1}{\pi \sqrt{(x - a)(b - x)}} 1873 * @f] 1874 * where @f$x >= a@f$ and @f$x <= b@f$. 1875 * 1876 * <table border=1 cellpadding=10 cellspacing=0> 1877 * <caption align=top>Distribution Statistics</caption> 1878 * <tr><td>Mean</td><td>@f$ (a + b) / 2 @f$</td></tr> 1879 * <tr><td>Variance</td><td>@f$ (b - a)^2 / 8 @f$</td></tr> 1880 * <tr><td>Range</td><td>@f$[a, b]@f$</td></tr> 1881 * </table> 1882 */ 1883 template<typename _RealType = double> 1884 class 1885 arcsine_distribution 1886 { 1887 static_assert(std::is_floating_point<_RealType>::value, 1888 "template argument not a floating point type"); 1889 1890 public: 1891 /** The type of the range of the distribution. */ 1892 typedef _RealType result_type; 1893 /** Parameter type. */ 1894 struct param_type 1895 { 1896 typedef arcsine_distribution<result_type> distribution_type; 1897 1898 param_type(result_type __a = result_type(0), 1899 result_type __b = result_type(1)) 1900 : _M_a(__a), _M_b(__b) 1901 { 1902 __glibcxx_assert(_M_a <= _M_b); 1903 } 1904 1905 result_type 1906 a() const 1907 { return _M_a; } 1908 1909 result_type 1910 b() const 1911 { return _M_b; } 1912 1913 friend bool 1914 operator==(const param_type& __p1, const param_type& __p2) 1915 { return __p1._M_a == __p2._M_a && __p1._M_b == __p2._M_b; } 1916 1917 private: 1918 void _M_initialize(); 1919 1920 result_type _M_a; 1921 result_type _M_b; 1922 }; 1923 1924 /** 1925 * @brief Constructors. 1926 */ 1927 explicit 1928 arcsine_distribution(result_type __a = result_type(0), 1929 result_type __b = result_type(1)) 1930 : _M_param(__a, __b), 1931 _M_ud(-1.5707963267948966192313216916397514L, 1932 +1.5707963267948966192313216916397514L) 1933 { } 1934 1935 explicit 1936 arcsine_distribution(const param_type& __p) 1937 : _M_param(__p), 1938 _M_ud(-1.5707963267948966192313216916397514L, 1939 +1.5707963267948966192313216916397514L) 1940 { } 1941 1942 /** 1943 * @brief Resets the distribution state. 1944 */ 1945 void 1946 reset() 1947 { _M_ud.reset(); } 1948 1949 /** 1950 * @brief Return the parameters of the distribution. 1951 */ 1952 result_type 1953 a() const 1954 { return _M_param.a(); } 1955 1956 result_type 1957 b() const 1958 { return _M_param.b(); } 1959 1960 /** 1961 * @brief Returns the parameter set of the distribution. 1962 */ 1963 param_type 1964 param() const 1965 { return _M_param; } 1966 1967 /** 1968 * @brief Sets the parameter set of the distribution. 1969 * @param __param The new parameter set of the distribution. 1970 */ 1971 void 1972 param(const param_type& __param) 1973 { _M_param = __param; } 1974 1975 /** 1976 * @brief Returns the greatest lower bound value of the distribution. 1977 */ 1978 result_type 1979 min() const 1980 { return this->a(); } 1981 1982 /** 1983 * @brief Returns the least upper bound value of the distribution. 1984 */ 1985 result_type 1986 max() const 1987 { return this->b(); } 1988 1989 /** 1990 * @brief Generating functions. 1991 */ 1992 template<typename _UniformRandomNumberGenerator> 1993 result_type 1994 operator()(_UniformRandomNumberGenerator& __urng) 1995 { 1996 result_type __x = std::sin(this->_M_ud(__urng)); 1997 return (__x * (this->b() - this->a()) 1998 + this->a() + this->b()) / result_type(2); 1999 } 2000 2001 template<typename _UniformRandomNumberGenerator> 2002 result_type 2003 operator()(_UniformRandomNumberGenerator& __urng, 2004 const param_type& __p) 2005 { 2006 result_type __x = std::sin(this->_M_ud(__urng)); 2007 return (__x * (__p.b() - __p.a()) 2008 + __p.a() + __p.b()) / result_type(2); 2009 } 2010 2011 template<typename _ForwardIterator, 2012 typename _UniformRandomNumberGenerator> 2013 void 2014 __generate(_ForwardIterator __f, _ForwardIterator __t, 2015 _UniformRandomNumberGenerator& __urng) 2016 { this->__generate(__f, __t, __urng, _M_param); } 2017 2018 template<typename _ForwardIterator, 2019 typename _UniformRandomNumberGenerator> 2020 void 2021 __generate(_ForwardIterator __f, _ForwardIterator __t, 2022 _UniformRandomNumberGenerator& __urng, 2023 const param_type& __p) 2024 { this->__generate_impl(__f, __t, __urng, __p); } 2025 2026 template<typename _UniformRandomNumberGenerator> 2027 void 2028 __generate(result_type* __f, result_type* __t, 2029 _UniformRandomNumberGenerator& __urng, 2030 const param_type& __p) 2031 { this->__generate_impl(__f, __t, __urng, __p); } 2032 2033 /** 2034 * @brief Return true if two arcsine distributions have 2035 * the same parameters and the sequences that would 2036 * be generated are equal. 2037 */ 2038 friend bool 2039 operator==(const arcsine_distribution& __d1, 2040 const arcsine_distribution& __d2) 2041 { return (__d1._M_param == __d2._M_param 2042 && __d1._M_ud == __d2._M_ud); } 2043 2044 /** 2045 * @brief Inserts a %arcsine_distribution random number distribution 2046 * @p __x into the output stream @p __os. 2047 * 2048 * @param __os An output stream. 2049 * @param __x A %arcsine_distribution random number distribution. 2050 * 2051 * @returns The output stream with the state of @p __x inserted or in 2052 * an error state. 2053 */ 2054 template<typename _RealType1, typename _CharT, typename _Traits> 2055 friend std::basic_ostream<_CharT, _Traits>& 2056 operator<<(std::basic_ostream<_CharT, _Traits>&, 2057 const arcsine_distribution<_RealType1>&); 2058 2059 /** 2060 * @brief Extracts a %arcsine_distribution random number distribution 2061 * @p __x from the input stream @p __is. 2062 * 2063 * @param __is An input stream. 2064 * @param __x A %arcsine_distribution random number 2065 * generator engine. 2066 * 2067 * @returns The input stream with @p __x extracted or in an error state. 2068 */ 2069 template<typename _RealType1, typename _CharT, typename _Traits> 2070 friend std::basic_istream<_CharT, _Traits>& 2071 operator>>(std::basic_istream<_CharT, _Traits>&, 2072 arcsine_distribution<_RealType1>&); 2073 2074 private: 2075 template<typename _ForwardIterator, 2076 typename _UniformRandomNumberGenerator> 2077 void 2078 __generate_impl(_ForwardIterator __f, _ForwardIterator __t, 2079 _UniformRandomNumberGenerator& __urng, 2080 const param_type& __p); 2081 2082 param_type _M_param; 2083 2084 std::uniform_real_distribution<result_type> _M_ud; 2085 }; 2086 2087 /** 2088 * @brief Return true if two arcsine distributions are not equal. 2089 */ 2090 template<typename _RealType> 2091 inline bool 2092 operator!=(const arcsine_distribution<_RealType>& __d1, 2093 const arcsine_distribution<_RealType>& __d2) 2094 { return !(__d1 == __d2); } 2095 2096 2097 /** 2098 * @brief A Hoyt continuous distribution for random numbers. 2099 * 2100 * The formula for the Hoyt probability density function is 2101 * @f[ 2102 * p(x|q,\omega) = \frac{(1 + q^2)x}{q\omega} 2103 * \exp\left(-\frac{(1 + q^2)^2 x^2}{4 q^2 \omega}\right) 2104 * I_0\left(\frac{(1 - q^4) x^2}{4 q^2 \omega}\right) 2105 * @f] 2106 * where @f$I_0(z)@f$ is the modified Bessel function of the first kind 2107 * of order 0 and @f$0 < q < 1@f$. 2108 * 2109 * <table border=1 cellpadding=10 cellspacing=0> 2110 * <caption align=top>Distribution Statistics</caption> 2111 * <tr><td>Mean</td><td>@f$ \sqrt{\frac{2}{\pi}} \sqrt{\frac{\omega}{1 + q^2}} 2112 * E(1 - q^2) @f$</td></tr> 2113 * <tr><td>Variance</td><td>@f$ \omega \left(1 - \frac{2E^2(1 - q^2)} 2114 * {\pi (1 + q^2)}\right) @f$</td></tr> 2115 * <tr><td>Range</td><td>@f$[0, \infty)@f$</td></tr> 2116 * </table> 2117 * where @f$E(x)@f$ is the elliptic function of the second kind. 2118 */ 2119 template<typename _RealType = double> 2120 class 2121 hoyt_distribution 2122 { 2123 static_assert(std::is_floating_point<_RealType>::value, 2124 "template argument not a floating point type"); 2125 2126 public: 2127 /** The type of the range of the distribution. */ 2128 typedef _RealType result_type; 2129 /** Parameter type. */ 2130 struct param_type 2131 { 2132 typedef hoyt_distribution<result_type> distribution_type; 2133 2134 param_type(result_type __q = result_type(0.5L), 2135 result_type __omega = result_type(1)) 2136 : _M_q(__q), _M_omega(__omega) 2137 { 2138 __glibcxx_assert(_M_q > result_type(0)); 2139 __glibcxx_assert(_M_q < result_type(1)); 2140 } 2141 2142 result_type 2143 q() const 2144 { return _M_q; } 2145 2146 result_type 2147 omega() const 2148 { return _M_omega; } 2149 2150 friend bool 2151 operator==(const param_type& __p1, const param_type& __p2) 2152 { return __p1._M_q == __p2._M_q 2153 && __p1._M_omega == __p2._M_omega; } 2154 2155 private: 2156 void _M_initialize(); 2157 2158 result_type _M_q; 2159 result_type _M_omega; 2160 }; 2161 2162 /** 2163 * @brief Constructors. 2164 */ 2165 explicit 2166 hoyt_distribution(result_type __q = result_type(0.5L), 2167 result_type __omega = result_type(1)) 2168 : _M_param(__q, __omega), 2169 _M_ad(result_type(0.5L) * (result_type(1) + __q * __q), 2170 result_type(0.5L) * (result_type(1) + __q * __q) 2171 / (__q * __q)), 2172 _M_ed(result_type(1)) 2173 { } 2174 2175 explicit 2176 hoyt_distribution(const param_type& __p) 2177 : _M_param(__p), 2178 _M_ad(result_type(0.5L) * (result_type(1) + __p.q() * __p.q()), 2179 result_type(0.5L) * (result_type(1) + __p.q() * __p.q()) 2180 / (__p.q() * __p.q())), 2181 _M_ed(result_type(1)) 2182 { } 2183 2184 /** 2185 * @brief Resets the distribution state. 2186 */ 2187 void 2188 reset() 2189 { 2190 _M_ad.reset(); 2191 _M_ed.reset(); 2192 } 2193 2194 /** 2195 * @brief Return the parameters of the distribution. 2196 */ 2197 result_type 2198 q() const 2199 { return _M_param.q(); } 2200 2201 result_type 2202 omega() const 2203 { return _M_param.omega(); } 2204 2205 /** 2206 * @brief Returns the parameter set of the distribution. 2207 */ 2208 param_type 2209 param() const 2210 { return _M_param; } 2211 2212 /** 2213 * @brief Sets the parameter set of the distribution. 2214 * @param __param The new parameter set of the distribution. 2215 */ 2216 void 2217 param(const param_type& __param) 2218 { _M_param = __param; } 2219 2220 /** 2221 * @brief Returns the greatest lower bound value of the distribution. 2222 */ 2223 result_type 2224 min() const 2225 { return result_type(0); } 2226 2227 /** 2228 * @brief Returns the least upper bound value of the distribution. 2229 */ 2230 result_type 2231 max() const 2232 { return std::numeric_limits<result_type>::max(); } 2233 2234 /** 2235 * @brief Generating functions. 2236 */ 2237 template<typename _UniformRandomNumberGenerator> 2238 result_type 2239 operator()(_UniformRandomNumberGenerator& __urng); 2240 2241 template<typename _UniformRandomNumberGenerator> 2242 result_type 2243 operator()(_UniformRandomNumberGenerator& __urng, 2244 const param_type& __p); 2245 2246 template<typename _ForwardIterator, 2247 typename _UniformRandomNumberGenerator> 2248 void 2249 __generate(_ForwardIterator __f, _ForwardIterator __t, 2250 _UniformRandomNumberGenerator& __urng) 2251 { this->__generate(__f, __t, __urng, _M_param); } 2252 2253 template<typename _ForwardIterator, 2254 typename _UniformRandomNumberGenerator> 2255 void 2256 __generate(_ForwardIterator __f, _ForwardIterator __t, 2257 _UniformRandomNumberGenerator& __urng, 2258 const param_type& __p) 2259 { this->__generate_impl(__f, __t, __urng, __p); } 2260 2261 template<typename _UniformRandomNumberGenerator> 2262 void 2263 __generate(result_type* __f, result_type* __t, 2264 _UniformRandomNumberGenerator& __urng, 2265 const param_type& __p) 2266 { this->__generate_impl(__f, __t, __urng, __p); } 2267 2268 /** 2269 * @brief Return true if two Hoyt distributions have 2270 * the same parameters and the sequences that would 2271 * be generated are equal. 2272 */ 2273 friend bool 2274 operator==(const hoyt_distribution& __d1, 2275 const hoyt_distribution& __d2) 2276 { return (__d1._M_param == __d2._M_param 2277 && __d1._M_ad == __d2._M_ad 2278 && __d1._M_ed == __d2._M_ed); } 2279 2280 /** 2281 * @brief Inserts a %hoyt_distribution random number distribution 2282 * @p __x into the output stream @p __os. 2283 * 2284 * @param __os An output stream. 2285 * @param __x A %hoyt_distribution random number distribution. 2286 * 2287 * @returns The output stream with the state of @p __x inserted or in 2288 * an error state. 2289 */ 2290 template<typename _RealType1, typename _CharT, typename _Traits> 2291 friend std::basic_ostream<_CharT, _Traits>& 2292 operator<<(std::basic_ostream<_CharT, _Traits>&, 2293 const hoyt_distribution<_RealType1>&); 2294 2295 /** 2296 * @brief Extracts a %hoyt_distribution random number distribution 2297 * @p __x from the input stream @p __is. 2298 * 2299 * @param __is An input stream. 2300 * @param __x A %hoyt_distribution random number 2301 * generator engine. 2302 * 2303 * @returns The input stream with @p __x extracted or in an error state. 2304 */ 2305 template<typename _RealType1, typename _CharT, typename _Traits> 2306 friend std::basic_istream<_CharT, _Traits>& 2307 operator>>(std::basic_istream<_CharT, _Traits>&, 2308 hoyt_distribution<_RealType1>&); 2309 2310 private: 2311 template<typename _ForwardIterator, 2312 typename _UniformRandomNumberGenerator> 2313 void 2314 __generate_impl(_ForwardIterator __f, _ForwardIterator __t, 2315 _UniformRandomNumberGenerator& __urng, 2316 const param_type& __p); 2317 2318 param_type _M_param; 2319 2320 __gnu_cxx::arcsine_distribution<result_type> _M_ad; 2321 std::exponential_distribution<result_type> _M_ed; 2322 }; 2323 2324 /** 2325 * @brief Return true if two Hoyt distributions are not equal. 2326 */ 2327 template<typename _RealType> 2328 inline bool 2329 operator!=(const hoyt_distribution<_RealType>& __d1, 2330 const hoyt_distribution<_RealType>& __d2) 2331 { return !(__d1 == __d2); } 2332 2333 2334 /** 2335 * @brief A triangular distribution for random numbers. 2336 * 2337 * The formula for the triangular probability density function is 2338 * @f[ 2339 * / 0 for x < a 2340 * p(x|a,b,c) = | \frac{2(x-a)}{(c-a)(b-a)} for a <= x <= b 2341 * | \frac{2(c-x)}{(c-a)(c-b)} for b < x <= c 2342 * \ 0 for c < x 2343 * @f] 2344 * 2345 * <table border=1 cellpadding=10 cellspacing=0> 2346 * <caption align=top>Distribution Statistics</caption> 2347 * <tr><td>Mean</td><td>@f$ \frac{a+b+c}{2} @f$</td></tr> 2348 * <tr><td>Variance</td><td>@f$ \frac{a^2+b^2+c^2-ab-ac-bc} 2349 * {18}@f$</td></tr> 2350 * <tr><td>Range</td><td>@f$[a, c]@f$</td></tr> 2351 * </table> 2352 */ 2353 template<typename _RealType = double> 2354 class triangular_distribution 2355 { 2356 static_assert(std::is_floating_point<_RealType>::value, 2357 "template argument not a floating point type"); 2358 2359 public: 2360 /** The type of the range of the distribution. */ 2361 typedef _RealType result_type; 2362 /** Parameter type. */ 2363 struct param_type 2364 { 2365 friend class triangular_distribution<_RealType>; 2366 2367 explicit 2368 param_type(_RealType __a = _RealType(0), 2369 _RealType __b = _RealType(0.5), 2370 _RealType __c = _RealType(1)) 2371 : _M_a(__a), _M_b(__b), _M_c(__c) 2372 { 2373 __glibcxx_assert(_M_a <= _M_b); 2374 __glibcxx_assert(_M_b <= _M_c); 2375 __glibcxx_assert(_M_a < _M_c); 2376 2377 _M_r_ab = (_M_b - _M_a) / (_M_c - _M_a); 2378 _M_f_ab_ac = (_M_b - _M_a) * (_M_c - _M_a); 2379 _M_f_bc_ac = (_M_c - _M_b) * (_M_c - _M_a); 2380 } 2381 2382 _RealType 2383 a() const 2384 { return _M_a; } 2385 2386 _RealType 2387 b() const 2388 { return _M_b; } 2389 2390 _RealType 2391 c() const 2392 { return _M_c; } 2393 2394 friend bool 2395 operator==(const param_type& __p1, const param_type& __p2) 2396 { return (__p1._M_a == __p2._M_a && __p1._M_b == __p2._M_b 2397 && __p1._M_c == __p2._M_c); } 2398 2399 private: 2400 2401 _RealType _M_a; 2402 _RealType _M_b; 2403 _RealType _M_c; 2404 _RealType _M_r_ab; 2405 _RealType _M_f_ab_ac; 2406 _RealType _M_f_bc_ac; 2407 }; 2408 2409 /** 2410 * @brief Constructs a triangle distribution with parameters 2411 * @f$ a @f$, @f$ b @f$ and @f$ c @f$. 2412 */ 2413 explicit 2414 triangular_distribution(result_type __a = result_type(0), 2415 result_type __b = result_type(0.5), 2416 result_type __c = result_type(1)) 2417 : _M_param(__a, __b, __c) 2418 { } 2419 2420 explicit 2421 triangular_distribution(const param_type& __p) 2422 : _M_param(__p) 2423 { } 2424 2425 /** 2426 * @brief Resets the distribution state. 2427 */ 2428 void 2429 reset() 2430 { } 2431 2432 /** 2433 * @brief Returns the @f$ a @f$ of the distribution. 2434 */ 2435 result_type 2436 a() const 2437 { return _M_param.a(); } 2438 2439 /** 2440 * @brief Returns the @f$ b @f$ of the distribution. 2441 */ 2442 result_type 2443 b() const 2444 { return _M_param.b(); } 2445 2446 /** 2447 * @brief Returns the @f$ c @f$ of the distribution. 2448 */ 2449 result_type 2450 c() const 2451 { return _M_param.c(); } 2452 2453 /** 2454 * @brief Returns the parameter set of the distribution. 2455 */ 2456 param_type 2457 param() const 2458 { return _M_param; } 2459 2460 /** 2461 * @brief Sets the parameter set of the distribution. 2462 * @param __param The new parameter set of the distribution. 2463 */ 2464 void 2465 param(const param_type& __param) 2466 { _M_param = __param; } 2467 2468 /** 2469 * @brief Returns the greatest lower bound value of the distribution. 2470 */ 2471 result_type 2472 min() const 2473 { return _M_param._M_a; } 2474 2475 /** 2476 * @brief Returns the least upper bound value of the distribution. 2477 */ 2478 result_type 2479 max() const 2480 { return _M_param._M_c; } 2481 2482 /** 2483 * @brief Generating functions. 2484 */ 2485 template<typename _UniformRandomNumberGenerator> 2486 result_type 2487 operator()(_UniformRandomNumberGenerator& __urng) 2488 { return this->operator()(__urng, _M_param); } 2489 2490 template<typename _UniformRandomNumberGenerator> 2491 result_type 2492 operator()(_UniformRandomNumberGenerator& __urng, 2493 const param_type& __p) 2494 { 2495 std::__detail::_Adaptor<_UniformRandomNumberGenerator, result_type> 2496 __aurng(__urng); 2497 result_type __rnd = __aurng(); 2498 if (__rnd <= __p._M_r_ab) 2499 return __p.a() + std::sqrt(__rnd * __p._M_f_ab_ac); 2500 else 2501 return __p.c() - std::sqrt((result_type(1) - __rnd) 2502 * __p._M_f_bc_ac); 2503 } 2504 2505 template<typename _ForwardIterator, 2506 typename _UniformRandomNumberGenerator> 2507 void 2508 __generate(_ForwardIterator __f, _ForwardIterator __t, 2509 _UniformRandomNumberGenerator& __urng) 2510 { this->__generate(__f, __t, __urng, _M_param); } 2511 2512 template<typename _ForwardIterator, 2513 typename _UniformRandomNumberGenerator> 2514 void 2515 __generate(_ForwardIterator __f, _ForwardIterator __t, 2516 _UniformRandomNumberGenerator& __urng, 2517 const param_type& __p) 2518 { this->__generate_impl(__f, __t, __urng, __p); } 2519 2520 template<typename _UniformRandomNumberGenerator> 2521 void 2522 __generate(result_type* __f, result_type* __t, 2523 _UniformRandomNumberGenerator& __urng, 2524 const param_type& __p) 2525 { this->__generate_impl(__f, __t, __urng, __p); } 2526 2527 /** 2528 * @brief Return true if two triangle distributions have the same 2529 * parameters and the sequences that would be generated 2530 * are equal. 2531 */ 2532 friend bool 2533 operator==(const triangular_distribution& __d1, 2534 const triangular_distribution& __d2) 2535 { return __d1._M_param == __d2._M_param; } 2536 2537 /** 2538 * @brief Inserts a %triangular_distribution random number distribution 2539 * @p __x into the output stream @p __os. 2540 * 2541 * @param __os An output stream. 2542 * @param __x A %triangular_distribution random number distribution. 2543 * 2544 * @returns The output stream with the state of @p __x inserted or in 2545 * an error state. 2546 */ 2547 template<typename _RealType1, typename _CharT, typename _Traits> 2548 friend std::basic_ostream<_CharT, _Traits>& 2549 operator<<(std::basic_ostream<_CharT, _Traits>& __os, 2550 const __gnu_cxx::triangular_distribution<_RealType1>& __x); 2551 2552 /** 2553 * @brief Extracts a %triangular_distribution random number distribution 2554 * @p __x from the input stream @p __is. 2555 * 2556 * @param __is An input stream. 2557 * @param __x A %triangular_distribution random number generator engine. 2558 * 2559 * @returns The input stream with @p __x extracted or in an error state. 2560 */ 2561 template<typename _RealType1, typename _CharT, typename _Traits> 2562 friend std::basic_istream<_CharT, _Traits>& 2563 operator>>(std::basic_istream<_CharT, _Traits>& __is, 2564 __gnu_cxx::triangular_distribution<_RealType1>& __x); 2565 2566 private: 2567 template<typename _ForwardIterator, 2568 typename _UniformRandomNumberGenerator> 2569 void 2570 __generate_impl(_ForwardIterator __f, _ForwardIterator __t, 2571 _UniformRandomNumberGenerator& __urng, 2572 const param_type& __p); 2573 2574 param_type _M_param; 2575 }; 2576 2577 /** 2578 * @brief Return true if two triangle distributions are different. 2579 */ 2580 template<typename _RealType> 2581 inline bool 2582 operator!=(const __gnu_cxx::triangular_distribution<_RealType>& __d1, 2583 const __gnu_cxx::triangular_distribution<_RealType>& __d2) 2584 { return !(__d1 == __d2); } 2585 2586 2587 /** 2588 * @brief A von Mises distribution for random numbers. 2589 * 2590 * The formula for the von Mises probability density function is 2591 * @f[ 2592 * p(x|\mu,\kappa) = \frac{e^{\kappa \cos(x-\mu)}} 2593 * {2\pi I_0(\kappa)} 2594 * @f] 2595 * 2596 * The generating functions use the method according to: 2597 * 2598 * D. J. Best and N. I. Fisher, 1979. "Efficient Simulation of the 2599 * von Mises Distribution", Journal of the Royal Statistical Society. 2600 * Series C (Applied Statistics), Vol. 28, No. 2, pp. 152-157. 2601 * 2602 * <table border=1 cellpadding=10 cellspacing=0> 2603 * <caption align=top>Distribution Statistics</caption> 2604 * <tr><td>Mean</td><td>@f$ \mu @f$</td></tr> 2605 * <tr><td>Variance</td><td>@f$ 1-I_1(\kappa)/I_0(\kappa) @f$</td></tr> 2606 * <tr><td>Range</td><td>@f$[-\pi, \pi]@f$</td></tr> 2607 * </table> 2608 */ 2609 template<typename _RealType = double> 2610 class von_mises_distribution 2611 { 2612 static_assert(std::is_floating_point<_RealType>::value, 2613 "template argument not a floating point type"); 2614 2615 public: 2616 /** The type of the range of the distribution. */ 2617 typedef _RealType result_type; 2618 /** Parameter type. */ 2619 struct param_type 2620 { 2621 friend class von_mises_distribution<_RealType>; 2622 2623 explicit 2624 param_type(_RealType __mu = _RealType(0), 2625 _RealType __kappa = _RealType(1)) 2626 : _M_mu(__mu), _M_kappa(__kappa) 2627 { 2628 const _RealType __pi = __gnu_cxx::__math_constants<_RealType>::__pi; 2629 __glibcxx_assert(_M_mu >= -__pi && _M_mu <= __pi); 2630 __glibcxx_assert(_M_kappa >= _RealType(0)); 2631 2632 auto __tau = std::sqrt(_RealType(4) * _M_kappa * _M_kappa 2633 + _RealType(1)) + _RealType(1); 2634 auto __rho = ((__tau - std::sqrt(_RealType(2) * __tau)) 2635 / (_RealType(2) * _M_kappa)); 2636 _M_r = (_RealType(1) + __rho * __rho) / (_RealType(2) * __rho); 2637 } 2638 2639 _RealType 2640 mu() const 2641 { return _M_mu; } 2642 2643 _RealType 2644 kappa() const 2645 { return _M_kappa; } 2646 2647 friend bool 2648 operator==(const param_type& __p1, const param_type& __p2) 2649 { return (__p1._M_mu == __p2._M_mu 2650 && __p1._M_kappa == __p2._M_kappa); } 2651 2652 private: 2653 _RealType _M_mu; 2654 _RealType _M_kappa; 2655 _RealType _M_r; 2656 }; 2657 2658 /** 2659 * @brief Constructs a von Mises distribution with parameters 2660 * @f$\mu@f$ and @f$\kappa@f$. 2661 */ 2662 explicit 2663 von_mises_distribution(result_type __mu = result_type(0), 2664 result_type __kappa = result_type(1)) 2665 : _M_param(__mu, __kappa) 2666 { } 2667 2668 explicit 2669 von_mises_distribution(const param_type& __p) 2670 : _M_param(__p) 2671 { } 2672 2673 /** 2674 * @brief Resets the distribution state. 2675 */ 2676 void 2677 reset() 2678 { } 2679 2680 /** 2681 * @brief Returns the @f$ \mu @f$ of the distribution. 2682 */ 2683 result_type 2684 mu() const 2685 { return _M_param.mu(); } 2686 2687 /** 2688 * @brief Returns the @f$ \kappa @f$ of the distribution. 2689 */ 2690 result_type 2691 kappa() const 2692 { return _M_param.kappa(); } 2693 2694 /** 2695 * @brief Returns the parameter set of the distribution. 2696 */ 2697 param_type 2698 param() const 2699 { return _M_param; } 2700 2701 /** 2702 * @brief Sets the parameter set of the distribution. 2703 * @param __param The new parameter set of the distribution. 2704 */ 2705 void 2706 param(const param_type& __param) 2707 { _M_param = __param; } 2708 2709 /** 2710 * @brief Returns the greatest lower bound value of the distribution. 2711 */ 2712 result_type 2713 min() const 2714 { 2715 return -__gnu_cxx::__math_constants<result_type>::__pi; 2716 } 2717 2718 /** 2719 * @brief Returns the least upper bound value of the distribution. 2720 */ 2721 result_type 2722 max() const 2723 { 2724 return __gnu_cxx::__math_constants<result_type>::__pi; 2725 } 2726 2727 /** 2728 * @brief Generating functions. 2729 */ 2730 template<typename _UniformRandomNumberGenerator> 2731 result_type 2732 operator()(_UniformRandomNumberGenerator& __urng) 2733 { return this->operator()(__urng, _M_param); } 2734 2735 template<typename _UniformRandomNumberGenerator> 2736 result_type 2737 operator()(_UniformRandomNumberGenerator& __urng, 2738 const param_type& __p); 2739 2740 template<typename _ForwardIterator, 2741 typename _UniformRandomNumberGenerator> 2742 void 2743 __generate(_ForwardIterator __f, _ForwardIterator __t, 2744 _UniformRandomNumberGenerator& __urng) 2745 { this->__generate(__f, __t, __urng, _M_param); } 2746 2747 template<typename _ForwardIterator, 2748 typename _UniformRandomNumberGenerator> 2749 void 2750 __generate(_ForwardIterator __f, _ForwardIterator __t, 2751 _UniformRandomNumberGenerator& __urng, 2752 const param_type& __p) 2753 { this->__generate_impl(__f, __t, __urng, __p); } 2754 2755 template<typename _UniformRandomNumberGenerator> 2756 void 2757 __generate(result_type* __f, result_type* __t, 2758 _UniformRandomNumberGenerator& __urng, 2759 const param_type& __p) 2760 { this->__generate_impl(__f, __t, __urng, __p); } 2761 2762 /** 2763 * @brief Return true if two von Mises distributions have the same 2764 * parameters and the sequences that would be generated 2765 * are equal. 2766 */ 2767 friend bool 2768 operator==(const von_mises_distribution& __d1, 2769 const von_mises_distribution& __d2) 2770 { return __d1._M_param == __d2._M_param; } 2771 2772 /** 2773 * @brief Inserts a %von_mises_distribution random number distribution 2774 * @p __x into the output stream @p __os. 2775 * 2776 * @param __os An output stream. 2777 * @param __x A %von_mises_distribution random number distribution. 2778 * 2779 * @returns The output stream with the state of @p __x inserted or in 2780 * an error state. 2781 */ 2782 template<typename _RealType1, typename _CharT, typename _Traits> 2783 friend std::basic_ostream<_CharT, _Traits>& 2784 operator<<(std::basic_ostream<_CharT, _Traits>& __os, 2785 const __gnu_cxx::von_mises_distribution<_RealType1>& __x); 2786 2787 /** 2788 * @brief Extracts a %von_mises_distribution random number distribution 2789 * @p __x from the input stream @p __is. 2790 * 2791 * @param __is An input stream. 2792 * @param __x A %von_mises_distribution random number generator engine. 2793 * 2794 * @returns The input stream with @p __x extracted or in an error state. 2795 */ 2796 template<typename _RealType1, typename _CharT, typename _Traits> 2797 friend std::basic_istream<_CharT, _Traits>& 2798 operator>>(std::basic_istream<_CharT, _Traits>& __is, 2799 __gnu_cxx::von_mises_distribution<_RealType1>& __x); 2800 2801 private: 2802 template<typename _ForwardIterator, 2803 typename _UniformRandomNumberGenerator> 2804 void 2805 __generate_impl(_ForwardIterator __f, _ForwardIterator __t, 2806 _UniformRandomNumberGenerator& __urng, 2807 const param_type& __p); 2808 2809 param_type _M_param; 2810 }; 2811 2812 /** 2813 * @brief Return true if two von Mises distributions are different. 2814 */ 2815 template<typename _RealType> 2816 inline bool 2817 operator!=(const __gnu_cxx::von_mises_distribution<_RealType>& __d1, 2818 const __gnu_cxx::von_mises_distribution<_RealType>& __d2) 2819 { return !(__d1 == __d2); } 2820 2821 2822 /** 2823 * @brief A discrete hypergeometric random number distribution. 2824 * 2825 * The hypergeometric distribution is a discrete probability distribution 2826 * that describes the probability of @p k successes in @p n draws @a without 2827 * replacement from a finite population of size @p N containing exactly @p K 2828 * successes. 2829 * 2830 * The formula for the hypergeometric probability density function is 2831 * @f[ 2832 * p(k|N,K,n) = \frac{\binom{K}{k} \binom{N-K}{n-k}}{\binom{N}{n}} 2833 * @f] 2834 * where @f$N@f$ is the total population of the distribution, 2835 * @f$K@f$ is the total population of the distribution. 2836 * 2837 * <table border=1 cellpadding=10 cellspacing=0> 2838 * <caption align=top>Distribution Statistics</caption> 2839 * <tr><td>Mean</td><td>@f$ n\frac{K}{N} @f$</td></tr> 2840 * <tr><td>Variance</td><td>@f$ n\frac{K}{N}\frac{N-K}{N}\frac{N-n}{N-1} 2841 * @f$</td></tr> 2842 * <tr><td>Range</td><td>@f$[max(0, n+K-N), min(K, n)]@f$</td></tr> 2843 * </table> 2844 */ 2845 template<typename _UIntType = unsigned int> 2846 class hypergeometric_distribution 2847 { 2848 static_assert(std::is_unsigned<_UIntType>::value, "template argument " 2849 "substituting _UIntType not an unsigned integral type"); 2850 2851 public: 2852 /** The type of the range of the distribution. */ 2853 typedef _UIntType result_type; 2854 2855 /** Parameter type. */ 2856 struct param_type 2857 { 2858 typedef hypergeometric_distribution<_UIntType> distribution_type; 2859 friend class hypergeometric_distribution<_UIntType>; 2860 2861 explicit 2862 param_type(result_type __N = 10, result_type __K = 5, 2863 result_type __n = 1) 2864 : _M_N{__N}, _M_K{__K}, _M_n{__n} 2865 { 2866 __glibcxx_assert(_M_N >= _M_K); 2867 __glibcxx_assert(_M_N >= _M_n); 2868 } 2869 2870 result_type 2871 total_size() const 2872 { return _M_N; } 2873 2874 result_type 2875 successful_size() const 2876 { return _M_K; } 2877 2878 result_type 2879 unsuccessful_size() const 2880 { return _M_N - _M_K; } 2881 2882 result_type 2883 total_draws() const 2884 { return _M_n; } 2885 2886 friend bool 2887 operator==(const param_type& __p1, const param_type& __p2) 2888 { return (__p1._M_N == __p2._M_N) 2889 && (__p1._M_K == __p2._M_K) 2890 && (__p1._M_n == __p2._M_n); } 2891 2892 private: 2893 2894 result_type _M_N; 2895 result_type _M_K; 2896 result_type _M_n; 2897 }; 2898 2899 // constructors and member function 2900 explicit 2901 hypergeometric_distribution(result_type __N = 10, result_type __K = 5, 2902 result_type __n = 1) 2903 : _M_param{__N, __K, __n} 2904 { } 2905 2906 explicit 2907 hypergeometric_distribution(const param_type& __p) 2908 : _M_param{__p} 2909 { } 2910 2911 /** 2912 * @brief Resets the distribution state. 2913 */ 2914 void 2915 reset() 2916 { } 2917 2918 /** 2919 * @brief Returns the distribution parameter @p N, 2920 * the total number of items. 2921 */ 2922 result_type 2923 total_size() const 2924 { return this->_M_param.total_size(); } 2925 2926 /** 2927 * @brief Returns the distribution parameter @p K, 2928 * the total number of successful items. 2929 */ 2930 result_type 2931 successful_size() const 2932 { return this->_M_param.successful_size(); } 2933 2934 /** 2935 * @brief Returns the total number of unsuccessful items @f$ N - K @f$. 2936 */ 2937 result_type 2938 unsuccessful_size() const 2939 { return this->_M_param.unsuccessful_size(); } 2940 2941 /** 2942 * @brief Returns the distribution parameter @p n, 2943 * the total number of draws. 2944 */ 2945 result_type 2946 total_draws() const 2947 { return this->_M_param.total_draws(); } 2948 2949 /** 2950 * @brief Returns the parameter set of the distribution. 2951 */ 2952 param_type 2953 param() const 2954 { return this->_M_param; } 2955 2956 /** 2957 * @brief Sets the parameter set of the distribution. 2958 * @param __param The new parameter set of the distribution. 2959 */ 2960 void 2961 param(const param_type& __param) 2962 { this->_M_param = __param; } 2963 2964 /** 2965 * @brief Returns the greatest lower bound value of the distribution. 2966 */ 2967 result_type 2968 min() const 2969 { 2970 using _IntType = typename std::make_signed<result_type>::type; 2971 return static_cast<result_type>(std::max(static_cast<_IntType>(0), 2972 static_cast<_IntType>(this->total_draws() 2973 - this->unsuccessful_size()))); 2974 } 2975 2976 /** 2977 * @brief Returns the least upper bound value of the distribution. 2978 */ 2979 result_type 2980 max() const 2981 { return std::min(this->successful_size(), this->total_draws()); } 2982 2983 /** 2984 * @brief Generating functions. 2985 */ 2986 template<typename _UniformRandomNumberGenerator> 2987 result_type 2988 operator()(_UniformRandomNumberGenerator& __urng) 2989 { return this->operator()(__urng, this->_M_param); } 2990 2991 template<typename _UniformRandomNumberGenerator> 2992 result_type 2993 operator()(_UniformRandomNumberGenerator& __urng, 2994 const param_type& __p); 2995 2996 template<typename _ForwardIterator, 2997 typename _UniformRandomNumberGenerator> 2998 void 2999 __generate(_ForwardIterator __f, _ForwardIterator __t, 3000 _UniformRandomNumberGenerator& __urng) 3001 { this->__generate(__f, __t, __urng, this->_M_param); } 3002 3003 template<typename _ForwardIterator, 3004 typename _UniformRandomNumberGenerator> 3005 void 3006 __generate(_ForwardIterator __f, _ForwardIterator __t, 3007 _UniformRandomNumberGenerator& __urng, 3008 const param_type& __p) 3009 { this->__generate_impl(__f, __t, __urng, __p); } 3010 3011 template<typename _UniformRandomNumberGenerator> 3012 void 3013 __generate(result_type* __f, result_type* __t, 3014 _UniformRandomNumberGenerator& __urng, 3015 const param_type& __p) 3016 { this->__generate_impl(__f, __t, __urng, __p); } 3017 3018 /** 3019 * @brief Return true if two hypergeometric distributions have the same 3020 * parameters and the sequences that would be generated 3021 * are equal. 3022 */ 3023 friend bool 3024 operator==(const hypergeometric_distribution& __d1, 3025 const hypergeometric_distribution& __d2) 3026 { return __d1._M_param == __d2._M_param; } 3027 3028 /** 3029 * @brief Inserts a %hypergeometric_distribution random number 3030 * distribution @p __x into the output stream @p __os. 3031 * 3032 * @param __os An output stream. 3033 * @param __x A %hypergeometric_distribution random number 3034 * distribution. 3035 * 3036 * @returns The output stream with the state of @p __x inserted or in 3037 * an error state. 3038 */ 3039 template<typename _UIntType1, typename _CharT, typename _Traits> 3040 friend std::basic_ostream<_CharT, _Traits>& 3041 operator<<(std::basic_ostream<_CharT, _Traits>& __os, 3042 const __gnu_cxx::hypergeometric_distribution<_UIntType1>& 3043 __x); 3044 3045 /** 3046 * @brief Extracts a %hypergeometric_distribution random number 3047 * distribution @p __x from the input stream @p __is. 3048 * 3049 * @param __is An input stream. 3050 * @param __x A %hypergeometric_distribution random number generator 3051 * distribution. 3052 * 3053 * @returns The input stream with @p __x extracted or in an error 3054 * state. 3055 */ 3056 template<typename _UIntType1, typename _CharT, typename _Traits> 3057 friend std::basic_istream<_CharT, _Traits>& 3058 operator>>(std::basic_istream<_CharT, _Traits>& __is, 3059 __gnu_cxx::hypergeometric_distribution<_UIntType1>& __x); 3060 3061 private: 3062 3063 template<typename _ForwardIterator, 3064 typename _UniformRandomNumberGenerator> 3065 void 3066 __generate_impl(_ForwardIterator __f, _ForwardIterator __t, 3067 _UniformRandomNumberGenerator& __urng, 3068 const param_type& __p); 3069 3070 param_type _M_param; 3071 }; 3072 3073 /** 3074 * @brief Return true if two hypergeometric distributions are different. 3075 */ 3076 template<typename _UIntType> 3077 inline bool 3078 operator!=(const __gnu_cxx::hypergeometric_distribution<_UIntType>& __d1, 3079 const __gnu_cxx::hypergeometric_distribution<_UIntType>& __d2) 3080 { return !(__d1 == __d2); } 3081 3082 /** 3083 * @brief A logistic continuous distribution for random numbers. 3084 * 3085 * The formula for the logistic probability density function is 3086 * @f[ 3087 * p(x|\a,\b) = \frac{e^{(x - a)/b}}{b[1 + e^{(x - a)/b}]^2} 3088 * @f] 3089 * where @f$b > 0@f$. 3090 * 3091 * The formula for the logistic probability function is 3092 * @f[ 3093 * cdf(x|\a,\b) = \frac{e^{(x - a)/b}}{1 + e^{(x - a)/b}} 3094 * @f] 3095 * where @f$b > 0@f$. 3096 * 3097 * <table border=1 cellpadding=10 cellspacing=0> 3098 * <caption align=top>Distribution Statistics</caption> 3099 * <tr><td>Mean</td><td>@f$a@f$</td></tr> 3100 * <tr><td>Variance</td><td>@f$b^2\pi^2/3@f$</td></tr> 3101 * <tr><td>Range</td><td>@f$[0, \infty)@f$</td></tr> 3102 * </table> 3103 */ 3104 template<typename _RealType = double> 3105 class 3106 logistic_distribution 3107 { 3108 static_assert(std::is_floating_point<_RealType>::value, 3109 "template argument not a floating point type"); 3110 3111 public: 3112 /** The type of the range of the distribution. */ 3113 typedef _RealType result_type; 3114 /** Parameter type. */ 3115 struct param_type 3116 { 3117 typedef logistic_distribution<result_type> distribution_type; 3118 3119 param_type(result_type __a = result_type(0), 3120 result_type __b = result_type(1)) 3121 : _M_a(__a), _M_b(__b) 3122 { 3123 __glibcxx_assert(_M_b > result_type(0)); 3124 } 3125 3126 result_type 3127 a() const 3128 { return _M_a; } 3129 3130 result_type 3131 b() const 3132 { return _M_b; } 3133 3134 friend bool 3135 operator==(const param_type& __p1, const param_type& __p2) 3136 { return __p1._M_a == __p2._M_a 3137 && __p1._M_b == __p2._M_b; } 3138 3139 private: 3140 void _M_initialize(); 3141 3142 result_type _M_a; 3143 result_type _M_b; 3144 }; 3145 3146 /** 3147 * @brief Constructors. 3148 */ 3149 explicit 3150 logistic_distribution(result_type __a = result_type(0), 3151 result_type __b = result_type(1)) 3152 : _M_param(__a, __b) 3153 { } 3154 3155 explicit 3156 logistic_distribution(const param_type& __p) 3157 : _M_param(__p) 3158 { } 3159 3160 /** 3161 * @brief Resets the distribution state. 3162 */ 3163 void 3164 reset() 3165 { } 3166 3167 /** 3168 * @brief Return the parameters of the distribution. 3169 */ 3170 result_type 3171 a() const 3172 { return _M_param.a(); } 3173 3174 result_type 3175 b() const 3176 { return _M_param.b(); } 3177 3178 /** 3179 * @brief Returns the parameter set of the distribution. 3180 */ 3181 param_type 3182 param() const 3183 { return _M_param; } 3184 3185 /** 3186 * @brief Sets the parameter set of the distribution. 3187 * @param __param The new parameter set of the distribution. 3188 */ 3189 void 3190 param(const param_type& __param) 3191 { _M_param = __param; } 3192 3193 /** 3194 * @brief Returns the greatest lower bound value of the distribution. 3195 */ 3196 result_type 3197 min() const 3198 { return -std::numeric_limits<result_type>::max(); } 3199 3200 /** 3201 * @brief Returns the least upper bound value of the distribution. 3202 */ 3203 result_type 3204 max() const 3205 { return std::numeric_limits<result_type>::max(); } 3206 3207 /** 3208 * @brief Generating functions. 3209 */ 3210 template<typename _UniformRandomNumberGenerator> 3211 result_type 3212 operator()(_UniformRandomNumberGenerator& __urng) 3213 { return this->operator()(__urng, this->_M_param); } 3214 3215 template<typename _UniformRandomNumberGenerator> 3216 result_type 3217 operator()(_UniformRandomNumberGenerator&, 3218 const param_type&); 3219 3220 template<typename _ForwardIterator, 3221 typename _UniformRandomNumberGenerator> 3222 void 3223 __generate(_ForwardIterator __f, _ForwardIterator __t, 3224 _UniformRandomNumberGenerator& __urng) 3225 { this->__generate(__f, __t, __urng, this->param()); } 3226 3227 template<typename _ForwardIterator, 3228 typename _UniformRandomNumberGenerator> 3229 void 3230 __generate(_ForwardIterator __f, _ForwardIterator __t, 3231 _UniformRandomNumberGenerator& __urng, 3232 const param_type& __p) 3233 { this->__generate_impl(__f, __t, __urng, __p); } 3234 3235 template<typename _UniformRandomNumberGenerator> 3236 void 3237 __generate(result_type* __f, result_type* __t, 3238 _UniformRandomNumberGenerator& __urng, 3239 const param_type& __p) 3240 { this->__generate_impl(__f, __t, __urng, __p); } 3241 3242 /** 3243 * @brief Return true if two logistic distributions have 3244 * the same parameters and the sequences that would 3245 * be generated are equal. 3246 */ 3247 template<typename _RealType1> 3248 friend bool 3249 operator==(const logistic_distribution<_RealType1>& __d1, 3250 const logistic_distribution<_RealType1>& __d2) 3251 { return __d1.param() == __d2.param(); } 3252 3253 /** 3254 * @brief Inserts a %logistic_distribution random number distribution 3255 * @p __x into the output stream @p __os. 3256 * 3257 * @param __os An output stream. 3258 * @param __x A %logistic_distribution random number distribution. 3259 * 3260 * @returns The output stream with the state of @p __x inserted or in 3261 * an error state. 3262 */ 3263 template<typename _RealType1, typename _CharT, typename _Traits> 3264 friend std::basic_ostream<_CharT, _Traits>& 3265 operator<<(std::basic_ostream<_CharT, _Traits>&, 3266 const logistic_distribution<_RealType1>&); 3267 3268 /** 3269 * @brief Extracts a %logistic_distribution random number distribution 3270 * @p __x from the input stream @p __is. 3271 * 3272 * @param __is An input stream. 3273 * @param __x A %logistic_distribution random number 3274 * generator engine. 3275 * 3276 * @returns The input stream with @p __x extracted or in an error state. 3277 */ 3278 template<typename _RealType1, typename _CharT, typename _Traits> 3279 friend std::basic_istream<_CharT, _Traits>& 3280 operator>>(std::basic_istream<_CharT, _Traits>&, 3281 logistic_distribution<_RealType1>&); 3282 3283 private: 3284 template<typename _ForwardIterator, 3285 typename _UniformRandomNumberGenerator> 3286 void 3287 __generate_impl(_ForwardIterator __f, _ForwardIterator __t, 3288 _UniformRandomNumberGenerator& __urng, 3289 const param_type& __p); 3290 3291 param_type _M_param; 3292 }; 3293 3294 /** 3295 * @brief Return true if two logistic distributions are not equal. 3296 */ 3297 template<typename _RealType1> 3298 inline bool 3299 operator!=(const logistic_distribution<_RealType1>& __d1, 3300 const logistic_distribution<_RealType1>& __d2) 3301 { return !(__d1 == __d2); } 3302 3303 3304 /** 3305 * @brief A distribution for random coordinates on a unit sphere. 3306 * 3307 * The method used in the generation function is attributed by Donald Knuth 3308 * to G. W. Brown, Modern Mathematics for the Engineer (1956). 3309 */ 3310 template<std::size_t _Dimen, typename _RealType = double> 3311 class uniform_on_sphere_distribution 3312 { 3313 static_assert(std::is_floating_point<_RealType>::value, 3314 "template argument not a floating point type"); 3315 static_assert(_Dimen != 0, "dimension is zero"); 3316 3317 public: 3318 /** The type of the range of the distribution. */ 3319 typedef std::array<_RealType, _Dimen> result_type; 3320 /** Parameter type. */ 3321 struct param_type 3322 { 3323 explicit 3324 param_type() 3325 { } 3326 3327 friend bool 3328 operator==(const param_type& __p1, const param_type& __p2) 3329 { return true; } 3330 }; 3331 3332 /** 3333 * @brief Constructs a uniform on sphere distribution. 3334 */ 3335 explicit 3336 uniform_on_sphere_distribution() 3337 : _M_param(), _M_nd() 3338 { } 3339 3340 explicit 3341 uniform_on_sphere_distribution(const param_type& __p) 3342 : _M_param(__p), _M_nd() 3343 { } 3344 3345 /** 3346 * @brief Resets the distribution state. 3347 */ 3348 void 3349 reset() 3350 { _M_nd.reset(); } 3351 3352 /** 3353 * @brief Returns the parameter set of the distribution. 3354 */ 3355 param_type 3356 param() const 3357 { return _M_param; } 3358 3359 /** 3360 * @brief Sets the parameter set of the distribution. 3361 * @param __param The new parameter set of the distribution. 3362 */ 3363 void 3364 param(const param_type& __param) 3365 { _M_param = __param; } 3366 3367 /** 3368 * @brief Returns the greatest lower bound value of the distribution. 3369 * This function makes no sense for this distribution. 3370 */ 3371 result_type 3372 min() const 3373 { 3374 result_type __res; 3375 __res.fill(0); 3376 return __res; 3377 } 3378 3379 /** 3380 * @brief Returns the least upper bound value of the distribution. 3381 * This function makes no sense for this distribution. 3382 */ 3383 result_type 3384 max() const 3385 { 3386 result_type __res; 3387 __res.fill(0); 3388 return __res; 3389 } 3390 3391 /** 3392 * @brief Generating functions. 3393 */ 3394 template<typename _UniformRandomNumberGenerator> 3395 result_type 3396 operator()(_UniformRandomNumberGenerator& __urng) 3397 { return this->operator()(__urng, _M_param); } 3398 3399 template<typename _UniformRandomNumberGenerator> 3400 result_type 3401 operator()(_UniformRandomNumberGenerator& __urng, 3402 const param_type& __p); 3403 3404 template<typename _ForwardIterator, 3405 typename _UniformRandomNumberGenerator> 3406 void 3407 __generate(_ForwardIterator __f, _ForwardIterator __t, 3408 _UniformRandomNumberGenerator& __urng) 3409 { this->__generate(__f, __t, __urng, this->param()); } 3410 3411 template<typename _ForwardIterator, 3412 typename _UniformRandomNumberGenerator> 3413 void 3414 __generate(_ForwardIterator __f, _ForwardIterator __t, 3415 _UniformRandomNumberGenerator& __urng, 3416 const param_type& __p) 3417 { this->__generate_impl(__f, __t, __urng, __p); } 3418 3419 template<typename _UniformRandomNumberGenerator> 3420 void 3421 __generate(result_type* __f, result_type* __t, 3422 _UniformRandomNumberGenerator& __urng, 3423 const param_type& __p) 3424 { this->__generate_impl(__f, __t, __urng, __p); } 3425 3426 /** 3427 * @brief Return true if two uniform on sphere distributions have 3428 * the same parameters and the sequences that would be 3429 * generated are equal. 3430 */ 3431 friend bool 3432 operator==(const uniform_on_sphere_distribution& __d1, 3433 const uniform_on_sphere_distribution& __d2) 3434 { return __d1._M_nd == __d2._M_nd; } 3435 3436 /** 3437 * @brief Inserts a %uniform_on_sphere_distribution random number 3438 * distribution @p __x into the output stream @p __os. 3439 * 3440 * @param __os An output stream. 3441 * @param __x A %uniform_on_sphere_distribution random number 3442 * distribution. 3443 * 3444 * @returns The output stream with the state of @p __x inserted or in 3445 * an error state. 3446 */ 3447 template<size_t _Dimen1, typename _RealType1, typename _CharT, 3448 typename _Traits> 3449 friend std::basic_ostream<_CharT, _Traits>& 3450 operator<<(std::basic_ostream<_CharT, _Traits>& __os, 3451 const __gnu_cxx::uniform_on_sphere_distribution<_Dimen1, 3452 _RealType1>& 3453 __x); 3454 3455 /** 3456 * @brief Extracts a %uniform_on_sphere_distribution random number 3457 * distribution 3458 * @p __x from the input stream @p __is. 3459 * 3460 * @param __is An input stream. 3461 * @param __x A %uniform_on_sphere_distribution random number 3462 * generator engine. 3463 * 3464 * @returns The input stream with @p __x extracted or in an error state. 3465 */ 3466 template<std::size_t _Dimen1, typename _RealType1, typename _CharT, 3467 typename _Traits> 3468 friend std::basic_istream<_CharT, _Traits>& 3469 operator>>(std::basic_istream<_CharT, _Traits>& __is, 3470 __gnu_cxx::uniform_on_sphere_distribution<_Dimen1, 3471 _RealType1>& __x); 3472 3473 private: 3474 template<typename _ForwardIterator, 3475 typename _UniformRandomNumberGenerator> 3476 void 3477 __generate_impl(_ForwardIterator __f, _ForwardIterator __t, 3478 _UniformRandomNumberGenerator& __urng, 3479 const param_type& __p); 3480 3481 param_type _M_param; 3482 std::normal_distribution<_RealType> _M_nd; 3483 }; 3484 3485 /** 3486 * @brief Return true if two uniform on sphere distributions are different. 3487 */ 3488 template<std::size_t _Dimen, typename _RealType> 3489 inline bool 3490 operator!=(const __gnu_cxx::uniform_on_sphere_distribution<_Dimen, 3491 _RealType>& __d1, 3492 const __gnu_cxx::uniform_on_sphere_distribution<_Dimen, 3493 _RealType>& __d2) 3494 { return !(__d1 == __d2); } 3495 3496_GLIBCXX_END_NAMESPACE_VERSION 3497} // namespace __gnu_cxx 3498 3499#include "ext/opt_random.h" 3500#include "random.tcc" 3501 3502#endif // _GLIBCXX_USE_C99_STDINT_TR1 && UINT32_C 3503 3504#endif // C++11 3505 3506#endif // _EXT_RANDOM 3507