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