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