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