| blob | eab3fb8ef3241d55c0a0df23aaa5e8a058a9b340 |
1 // Random number extensions -*- C++ -*-
3 // Copyright (C) 2012-2014 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.
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.
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.
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/>.
25 /** @file ext/random
26 * This file is a GNU extension to the Standard C++ Library.
27 */
29 #ifndef _EXT_RANDOM
30 #define _EXT_RANDOM 1
32 #pragma GCC system_header
34 #if __cplusplus < 201103L
35 # include <bits/c++0x_warning.h>
36 #else
38 #include <random>
39 #include <array>
40 #include <ext/cmath>
41 #ifdef __SSE2__
42 # include <x86intrin.h>
43 #endif
45 #ifdef _GLIBCXX_USE_C99_STDINT_TR1
47 namespace __gnu_cxx _GLIBCXX_VISIBILITY(default)
48 {
49 _GLIBCXX_BEGIN_NAMESPACE_VERSION
51 #if __BYTE_ORDER__ == __ORDER_LITTLE_ENDIAN__
53 /* Mersenne twister implementation optimized for vector operations.
54 *
55 * Reference: http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/SFMT/
56 */
57 template<typename _UIntType, size_t __m,
58 size_t __pos1, size_t __sl1, size_t __sl2,
59 size_t __sr1, size_t __sr2,
60 uint32_t __msk1, uint32_t __msk2,
61 uint32_t __msk3, uint32_t __msk4,
62 uint32_t __parity1, uint32_t __parity2,
63 uint32_t __parity3, uint32_t __parity4>
64 class simd_fast_mersenne_twister_engine
65 {
66 static_assert(std::is_unsigned<_UIntType>::value, "template argument "
67 "substituting _UIntType not an unsigned integral type");
68 static_assert(__sr1 < 32, "first right shift too large");
69 static_assert(__sr2 < 16, "second right shift too large");
70 static_assert(__sl1 < 32, "first left shift too large");
71 static_assert(__sl2 < 16, "second left shift too large");
73 public:
74 typedef _UIntType result_type;
76 private:
77 static constexpr size_t m_w = sizeof(result_type) * 8;
78 static constexpr size_t _M_nstate = __m / 128 + 1;
79 static constexpr size_t _M_nstate32 = _M_nstate * 4;
81 static_assert(std::is_unsigned<_UIntType>::value, "template argument "
82 "substituting _UIntType not an unsigned integral type");
83 static_assert(__pos1 < _M_nstate, "POS1 not smaller than state size");
84 static_assert(16 % sizeof(_UIntType) == 0,
85 "UIntType size must divide 16");
87 public:
88 static constexpr size_t state_size = _M_nstate * (16
89 / sizeof(result_type));
90 static constexpr result_type default_seed = 5489u;
92 // constructors and member function
93 explicit
94 simd_fast_mersenne_twister_engine(result_type __sd = default_seed)
95 { seed(__sd); }
97 template<typename _Sseq, typename = typename
98 std::enable_if<!std::is_same<_Sseq,
99 simd_fast_mersenne_twister_engine>::value>
100 ::type>
101 explicit
102 simd_fast_mersenne_twister_engine(_Sseq& __q)
103 { seed(__q); }
105 void
106 seed(result_type __sd = default_seed);
108 template<typename _Sseq>
109 typename std::enable_if<std::is_class<_Sseq>::value>::type
110 seed(_Sseq& __q);
112 static constexpr result_type
113 min()
114 { return 0; };
116 static constexpr result_type
117 max()
118 { return std::numeric_limits<result_type>::max(); }
120 void
121 discard(unsigned long long __z);
123 result_type
124 operator()()
125 {
126 if (__builtin_expect(_M_pos >= state_size, 0))
127 _M_gen_rand();
129 return _M_stateT[_M_pos++];
130 }
132 template<typename _UIntType_2, size_t __m_2,
133 size_t __pos1_2, size_t __sl1_2, size_t __sl2_2,
134 size_t __sr1_2, size_t __sr2_2,
135 uint32_t __msk1_2, uint32_t __msk2_2,
136 uint32_t __msk3_2, uint32_t __msk4_2,
137 uint32_t __parity1_2, uint32_t __parity2_2,
138 uint32_t __parity3_2, uint32_t __parity4_2>
139 friend bool
140 operator==(const simd_fast_mersenne_twister_engine<_UIntType_2,
141 __m_2, __pos1_2, __sl1_2, __sl2_2, __sr1_2, __sr2_2,
142 __msk1_2, __msk2_2, __msk3_2, __msk4_2,
143 __parity1_2, __parity2_2, __parity3_2, __parity4_2>& __lhs,
144 const simd_fast_mersenne_twister_engine<_UIntType_2,
145 __m_2, __pos1_2, __sl1_2, __sl2_2, __sr1_2, __sr2_2,
146 __msk1_2, __msk2_2, __msk3_2, __msk4_2,
147 __parity1_2, __parity2_2, __parity3_2, __parity4_2>& __rhs);
149 template<typename _UIntType_2, size_t __m_2,
150 size_t __pos1_2, size_t __sl1_2, size_t __sl2_2,
151 size_t __sr1_2, size_t __sr2_2,
152 uint32_t __msk1_2, uint32_t __msk2_2,
153 uint32_t __msk3_2, uint32_t __msk4_2,
154 uint32_t __parity1_2, uint32_t __parity2_2,
155 uint32_t __parity3_2, uint32_t __parity4_2,
156 typename _CharT, typename _Traits>
157 friend std::basic_ostream<_CharT, _Traits>&
158 operator<<(std::basic_ostream<_CharT, _Traits>& __os,
159 const __gnu_cxx::simd_fast_mersenne_twister_engine
160 <_UIntType_2,
161 __m_2, __pos1_2, __sl1_2, __sl2_2, __sr1_2, __sr2_2,
162 __msk1_2, __msk2_2, __msk3_2, __msk4_2,
163 __parity1_2, __parity2_2, __parity3_2, __parity4_2>& __x);
165 template<typename _UIntType_2, size_t __m_2,
166 size_t __pos1_2, size_t __sl1_2, size_t __sl2_2,
167 size_t __sr1_2, size_t __sr2_2,
168 uint32_t __msk1_2, uint32_t __msk2_2,
169 uint32_t __msk3_2, uint32_t __msk4_2,
170 uint32_t __parity1_2, uint32_t __parity2_2,
171 uint32_t __parity3_2, uint32_t __parity4_2,
172 typename _CharT, typename _Traits>
173 friend std::basic_istream<_CharT, _Traits>&
174 operator>>(std::basic_istream<_CharT, _Traits>& __is,
175 __gnu_cxx::simd_fast_mersenne_twister_engine<_UIntType_2,
176 __m_2, __pos1_2, __sl1_2, __sl2_2, __sr1_2, __sr2_2,
177 __msk1_2, __msk2_2, __msk3_2, __msk4_2,
178 __parity1_2, __parity2_2, __parity3_2, __parity4_2>& __x);
180 private:
181 union
182 {
183 #ifdef __SSE2__
184 __m128i _M_state[_M_nstate];
185 #endif
186 uint32_t _M_state32[_M_nstate32];
187 result_type _M_stateT[state_size];
188 } __attribute__ ((__aligned__ (16)));
189 size_t _M_pos;
191 void _M_gen_rand(void);
192 void _M_period_certification();
193 };
196 template<typename _UIntType, size_t __m,
197 size_t __pos1, size_t __sl1, size_t __sl2,
198 size_t __sr1, size_t __sr2,
199 uint32_t __msk1, uint32_t __msk2,
200 uint32_t __msk3, uint32_t __msk4,
201 uint32_t __parity1, uint32_t __parity2,
202 uint32_t __parity3, uint32_t __parity4>
203 inline bool
204 operator!=(const __gnu_cxx::simd_fast_mersenne_twister_engine<_UIntType,
205 __m, __pos1, __sl1, __sl2, __sr1, __sr2, __msk1, __msk2, __msk3,
206 __msk4, __parity1, __parity2, __parity3, __parity4>& __lhs,
207 const __gnu_cxx::simd_fast_mersenne_twister_engine<_UIntType,
208 __m, __pos1, __sl1, __sl2, __sr1, __sr2, __msk1, __msk2, __msk3,
209 __msk4, __parity1, __parity2, __parity3, __parity4>& __rhs)
210 { return !(__lhs == __rhs); }
213 /* Definitions for the SIMD-oriented Fast Mersenne Twister as defined
214 * in the C implementation by Daito and Matsumoto, as both a 32-bit
215 * and 64-bit version.
216 */
217 typedef simd_fast_mersenne_twister_engine<uint32_t, 607, 2,
218 15, 3, 13, 3,
219 0xfdff37ffU, 0xef7f3f7dU,
220 0xff777b7dU, 0x7ff7fb2fU,
221 0x00000001U, 0x00000000U,
222 0x00000000U, 0x5986f054U>
223 sfmt607;
225 typedef simd_fast_mersenne_twister_engine<uint64_t, 607, 2,
226 15, 3, 13, 3,
227 0xfdff37ffU, 0xef7f3f7dU,
228 0xff777b7dU, 0x7ff7fb2fU,
229 0x00000001U, 0x00000000U,
230 0x00000000U, 0x5986f054U>
231 sfmt607_64;
234 typedef simd_fast_mersenne_twister_engine<uint32_t, 1279, 7,
235 14, 3, 5, 1,
236 0xf7fefffdU, 0x7fefcfffU,
237 0xaff3ef3fU, 0xb5ffff7fU,
238 0x00000001U, 0x00000000U,
239 0x00000000U, 0x20000000U>
240 sfmt1279;
242 typedef simd_fast_mersenne_twister_engine<uint64_t, 1279, 7,
243 14, 3, 5, 1,
244 0xf7fefffdU, 0x7fefcfffU,
245 0xaff3ef3fU, 0xb5ffff7fU,
246 0x00000001U, 0x00000000U,
247 0x00000000U, 0x20000000U>
248 sfmt1279_64;
251 typedef simd_fast_mersenne_twister_engine<uint32_t, 2281, 12,
252 19, 1, 5, 1,
253 0xbff7ffbfU, 0xfdfffffeU,
254 0xf7ffef7fU, 0xf2f7cbbfU,
255 0x00000001U, 0x00000000U,
256 0x00000000U, 0x41dfa600U>
257 sfmt2281;
259 typedef simd_fast_mersenne_twister_engine<uint64_t, 2281, 12,
260 19, 1, 5, 1,
261 0xbff7ffbfU, 0xfdfffffeU,
262 0xf7ffef7fU, 0xf2f7cbbfU,
263 0x00000001U, 0x00000000U,
264 0x00000000U, 0x41dfa600U>
265 sfmt2281_64;
268 typedef simd_fast_mersenne_twister_engine<uint32_t, 4253, 17,
269 20, 1, 7, 1,
270 0x9f7bffffU, 0x9fffff5fU,
271 0x3efffffbU, 0xfffff7bbU,
272 0xa8000001U, 0xaf5390a3U,
273 0xb740b3f8U, 0x6c11486dU>
274 sfmt4253;
276 typedef simd_fast_mersenne_twister_engine<uint64_t, 4253, 17,
277 20, 1, 7, 1,
278 0x9f7bffffU, 0x9fffff5fU,
279 0x3efffffbU, 0xfffff7bbU,
280 0xa8000001U, 0xaf5390a3U,
281 0xb740b3f8U, 0x6c11486dU>
282 sfmt4253_64;
285 typedef simd_fast_mersenne_twister_engine<uint32_t, 11213, 68,
286 14, 3, 7, 3,
287 0xeffff7fbU, 0xffffffefU,
288 0xdfdfbfffU, 0x7fffdbfdU,
289 0x00000001U, 0x00000000U,
290 0xe8148000U, 0xd0c7afa3U>
291 sfmt11213;
293 typedef simd_fast_mersenne_twister_engine<uint64_t, 11213, 68,
294 14, 3, 7, 3,
295 0xeffff7fbU, 0xffffffefU,
296 0xdfdfbfffU, 0x7fffdbfdU,
297 0x00000001U, 0x00000000U,
298 0xe8148000U, 0xd0c7afa3U>
299 sfmt11213_64;
302 typedef simd_fast_mersenne_twister_engine<uint32_t, 19937, 122,
303 18, 1, 11, 1,
304 0xdfffffefU, 0xddfecb7fU,
305 0xbffaffffU, 0xbffffff6U,
306 0x00000001U, 0x00000000U,
307 0x00000000U, 0x13c9e684U>
308 sfmt19937;
310 typedef simd_fast_mersenne_twister_engine<uint64_t, 19937, 122,
311 18, 1, 11, 1,
312 0xdfffffefU, 0xddfecb7fU,
313 0xbffaffffU, 0xbffffff6U,
314 0x00000001U, 0x00000000U,
315 0x00000000U, 0x13c9e684U>
316 sfmt19937_64;
319 typedef simd_fast_mersenne_twister_engine<uint32_t, 44497, 330,
320 5, 3, 9, 3,
321 0xeffffffbU, 0xdfbebfffU,
322 0xbfbf7befU, 0x9ffd7bffU,
323 0x00000001U, 0x00000000U,
324 0xa3ac4000U, 0xecc1327aU>
325 sfmt44497;
327 typedef simd_fast_mersenne_twister_engine<uint64_t, 44497, 330,
328 5, 3, 9, 3,
329 0xeffffffbU, 0xdfbebfffU,
330 0xbfbf7befU, 0x9ffd7bffU,
331 0x00000001U, 0x00000000U,
332 0xa3ac4000U, 0xecc1327aU>
333 sfmt44497_64;
336 typedef simd_fast_mersenne_twister_engine<uint32_t, 86243, 366,
337 6, 7, 19, 1,
338 0xfdbffbffU, 0xbff7ff3fU,
339 0xfd77efffU, 0xbf9ff3ffU,
340 0x00000001U, 0x00000000U,
341 0x00000000U, 0xe9528d85U>
342 sfmt86243;
344 typedef simd_fast_mersenne_twister_engine<uint64_t, 86243, 366,
345 6, 7, 19, 1,
346 0xfdbffbffU, 0xbff7ff3fU,
347 0xfd77efffU, 0xbf9ff3ffU,
348 0x00000001U, 0x00000000U,
349 0x00000000U, 0xe9528d85U>
350 sfmt86243_64;
353 typedef simd_fast_mersenne_twister_engine<uint32_t, 132049, 110,
354 19, 1, 21, 1,
355 0xffffbb5fU, 0xfb6ebf95U,
356 0xfffefffaU, 0xcff77fffU,
357 0x00000001U, 0x00000000U,
358 0xcb520000U, 0xc7e91c7dU>
359 sfmt132049;
361 typedef simd_fast_mersenne_twister_engine<uint64_t, 132049, 110,
362 19, 1, 21, 1,
363 0xffffbb5fU, 0xfb6ebf95U,
364 0xfffefffaU, 0xcff77fffU,
365 0x00000001U, 0x00000000U,
366 0xcb520000U, 0xc7e91c7dU>
367 sfmt132049_64;
370 typedef simd_fast_mersenne_twister_engine<uint32_t, 216091, 627,
371 11, 3, 10, 1,
372 0xbff7bff7U, 0xbfffffffU,
373 0xbffffa7fU, 0xffddfbfbU,
374 0xf8000001U, 0x89e80709U,
375 0x3bd2b64bU, 0x0c64b1e4U>
376 sfmt216091;
378 typedef simd_fast_mersenne_twister_engine<uint64_t, 216091, 627,
379 11, 3, 10, 1,
380 0xbff7bff7U, 0xbfffffffU,
381 0xbffffa7fU, 0xffddfbfbU,
382 0xf8000001U, 0x89e80709U,
383 0x3bd2b64bU, 0x0c64b1e4U>
384 sfmt216091_64;
386 #endif // __BYTE_ORDER__ == __ORDER_LITTLE_ENDIAN__
388 /**
389 * @brief A beta continuous distribution for random numbers.
390 *
391 * The formula for the beta probability density function is:
392 * @f[
393 * p(x|\alpha,\beta) = \frac{1}{B(\alpha,\beta)}
394 * x^{\alpha - 1} (1 - x)^{\beta - 1}
395 * @f]
396 */
397 template<typename _RealType = double>
398 class beta_distribution
399 {
400 static_assert(std::is_floating_point<_RealType>::value,
401 "template argument not a floating point type");
403 public:
404 /** The type of the range of the distribution. */
405 typedef _RealType result_type;
406 /** Parameter type. */
407 struct param_type
408 {
409 typedef beta_distribution<_RealType> distribution_type;
410 friend class beta_distribution<_RealType>;
412 explicit
413 param_type(_RealType __alpha_val = _RealType(1),
414 _RealType __beta_val = _RealType(1))
415 : _M_alpha(__alpha_val), _M_beta(__beta_val)
416 {
417 _GLIBCXX_DEBUG_ASSERT(_M_alpha > _RealType(0));
418 _GLIBCXX_DEBUG_ASSERT(_M_beta > _RealType(0));
419 }
421 _RealType
422 alpha() const
423 { return _M_alpha; }
425 _RealType
426 beta() const
427 { return _M_beta; }
429 friend bool
430 operator==(const param_type& __p1, const param_type& __p2)
431 { return (__p1._M_alpha == __p2._M_alpha
432 && __p1._M_beta == __p2._M_beta); }
434 private:
435 void
436 _M_initialize();
438 _RealType _M_alpha;
439 _RealType _M_beta;
440 };
442 public:
443 /**
444 * @brief Constructs a beta distribution with parameters
445 * @f$\alpha@f$ and @f$\beta@f$.
446 */
447 explicit
448 beta_distribution(_RealType __alpha_val = _RealType(1),
449 _RealType __beta_val = _RealType(1))
450 : _M_param(__alpha_val, __beta_val)
451 { }
453 explicit
454 beta_distribution(const param_type& __p)
455 : _M_param(__p)
456 { }
458 /**
459 * @brief Resets the distribution state.
460 */
461 void
462 reset()
463 { }
465 /**
466 * @brief Returns the @f$\alpha@f$ of the distribution.
467 */
468 _RealType
469 alpha() const
470 { return _M_param.alpha(); }
472 /**
473 * @brief Returns the @f$\beta@f$ of the distribution.
474 */
475 _RealType
476 beta() const
477 { return _M_param.beta(); }
479 /**
480 * @brief Returns the parameter set of the distribution.
481 */
482 param_type
483 param() const
484 { return _M_param; }
486 /**
487 * @brief Sets the parameter set of the distribution.
488 * @param __param The new parameter set of the distribution.
489 */
490 void
491 param(const param_type& __param)
492 { _M_param = __param; }
494 /**
495 * @brief Returns the greatest lower bound value of the distribution.
496 */
497 result_type
498 min() const
499 { return result_type(0); }
501 /**
502 * @brief Returns the least upper bound value of the distribution.
503 */
504 result_type
505 max() const
506 { return result_type(1); }
508 /**
509 * @brief Generating functions.
510 */
511 template<typename _UniformRandomNumberGenerator>
512 result_type
513 operator()(_UniformRandomNumberGenerator& __urng)
514 { return this->operator()(__urng, _M_param); }
516 template<typename _UniformRandomNumberGenerator>
517 result_type
518 operator()(_UniformRandomNumberGenerator& __urng,
519 const param_type& __p);
521 template<typename _ForwardIterator,
522 typename _UniformRandomNumberGenerator>
523 void
524 __generate(_ForwardIterator __f, _ForwardIterator __t,
525 _UniformRandomNumberGenerator& __urng)
526 { this->__generate(__f, __t, __urng, _M_param); }
528 template<typename _ForwardIterator,
529 typename _UniformRandomNumberGenerator>
530 void
531 __generate(_ForwardIterator __f, _ForwardIterator __t,
532 _UniformRandomNumberGenerator& __urng,
533 const param_type& __p)
534 { this->__generate_impl(__f, __t, __urng, __p); }
536 template<typename _UniformRandomNumberGenerator>
537 void
538 __generate(result_type* __f, result_type* __t,
539 _UniformRandomNumberGenerator& __urng,
540 const param_type& __p)
541 { this->__generate_impl(__f, __t, __urng, __p); }
543 /**
544 * @brief Return true if two beta distributions have the same
545 * parameters and the sequences that would be generated
546 * are equal.
547 */
548 friend bool
549 operator==(const beta_distribution& __d1,
550 const beta_distribution& __d2)
551 { return __d1._M_param == __d2._M_param; }
553 /**
554 * @brief Inserts a %beta_distribution random number distribution
555 * @p __x into the output stream @p __os.
556 *
557 * @param __os An output stream.
558 * @param __x A %beta_distribution random number distribution.
559 *
560 * @returns The output stream with the state of @p __x inserted or in
561 * an error state.
562 */
563 template<typename _RealType1, typename _CharT, typename _Traits>
564 friend std::basic_ostream<_CharT, _Traits>&
565 operator<<(std::basic_ostream<_CharT, _Traits>& __os,
566 const __gnu_cxx::beta_distribution<_RealType1>& __x);
568 /**
569 * @brief Extracts a %beta_distribution random number distribution
570 * @p __x from the input stream @p __is.
571 *
572 * @param __is An input stream.
573 * @param __x A %beta_distribution random number generator engine.
574 *
575 * @returns The input stream with @p __x extracted or in an error state.
576 */
577 template<typename _RealType1, typename _CharT, typename _Traits>
578 friend std::basic_istream<_CharT, _Traits>&
579 operator>>(std::basic_istream<_CharT, _Traits>& __is,
580 __gnu_cxx::beta_distribution<_RealType1>& __x);
582 private:
583 template<typename _ForwardIterator,
584 typename _UniformRandomNumberGenerator>
585 void
586 __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
587 _UniformRandomNumberGenerator& __urng,
588 const param_type& __p);
590 param_type _M_param;
591 };
593 /**
594 * @brief Return true if two beta distributions are different.
595 */
596 template<typename _RealType>
597 inline bool
598 operator!=(const __gnu_cxx::beta_distribution<_RealType>& __d1,
599 const __gnu_cxx::beta_distribution<_RealType>& __d2)
600 { return !(__d1 == __d2); }
603 /**
604 * @brief A multi-variate normal continuous distribution for random numbers.
605 *
606 * The formula for the normal probability density function is
607 * @f[
608 * p(\overrightarrow{x}|\overrightarrow{\mu },\Sigma) =
609 * \frac{1}{\sqrt{(2\pi )^k\det(\Sigma))}}
610 * e^{-\frac{1}{2}(\overrightarrow{x}-\overrightarrow{\mu})^\text{T}
611 * \Sigma ^{-1}(\overrightarrow{x}-\overrightarrow{\mu})}
612 * @f]
613 *
614 * where @f$\overrightarrow{x}@f$ and @f$\overrightarrow{\mu}@f$ are
615 * vectors of dimension @f$k@f$ and @f$\Sigma@f$ is the covariance
616 * matrix (which must be positive-definite).
617 */
618 template<std::size_t _Dimen, typename _RealType = double>
619 class normal_mv_distribution
620 {
621 static_assert(std::is_floating_point<_RealType>::value,
622 "template argument not a floating point type");
623 static_assert(_Dimen != 0, "dimension is zero");
625 public:
626 /** The type of the range of the distribution. */
627 typedef std::array<_RealType, _Dimen> result_type;
628 /** Parameter type. */
629 class param_type
630 {
631 static constexpr size_t _M_t_size = _Dimen * (_Dimen + 1) / 2;
633 public:
634 typedef normal_mv_distribution<_Dimen, _RealType> distribution_type;
635 friend class normal_mv_distribution<_Dimen, _RealType>;
637 param_type()
638 {
639 std::fill(_M_mean.begin(), _M_mean.end(), _RealType(0));
640 auto __it = _M_t.begin();
641 for (size_t __i = 0; __i < _Dimen; ++__i)
642 {
643 std::fill_n(__it, __i, _RealType(0));
644 __it += __i;
645 *__it++ = _RealType(1);
646 }
647 }
649 template<typename _ForwardIterator1, typename _ForwardIterator2>
650 param_type(_ForwardIterator1 __meanbegin,
651 _ForwardIterator1 __meanend,
652 _ForwardIterator2 __varcovbegin,
653 _ForwardIterator2 __varcovend)
654 {
655 __glibcxx_function_requires(_ForwardIteratorConcept<
656 _ForwardIterator1>)
657 __glibcxx_function_requires(_ForwardIteratorConcept<
658 _ForwardIterator2>)
659 _GLIBCXX_DEBUG_ASSERT(std::distance(__meanbegin, __meanend)
660 <= _Dimen);
661 const auto __dist = std::distance(__varcovbegin, __varcovend);
662 _GLIBCXX_DEBUG_ASSERT(__dist == _Dimen * _Dimen
663 || __dist == _Dimen * (_Dimen + 1) / 2
664 || __dist == _Dimen);
666 if (__dist == _Dimen * _Dimen)
667 _M_init_full(__meanbegin, __meanend, __varcovbegin, __varcovend);
668 else if (__dist == _Dimen * (_Dimen + 1) / 2)
669 _M_init_lower(__meanbegin, __meanend, __varcovbegin, __varcovend);
670 else
671 _M_init_diagonal(__meanbegin, __meanend,
672 __varcovbegin, __varcovend);
673 }
675 param_type(std::initializer_list<_RealType> __mean,
676 std::initializer_list<_RealType> __varcov)
677 {
678 _GLIBCXX_DEBUG_ASSERT(__mean.size() <= _Dimen);
679 _GLIBCXX_DEBUG_ASSERT(__varcov.size() == _Dimen * _Dimen
680 || __varcov.size() == _Dimen * (_Dimen + 1) / 2
681 || __varcov.size() == _Dimen);
683 if (__varcov.size() == _Dimen * _Dimen)
684 _M_init_full(__mean.begin(), __mean.end(),
685 __varcov.begin(), __varcov.end());
686 else if (__varcov.size() == _Dimen * (_Dimen + 1) / 2)
687 _M_init_lower(__mean.begin(), __mean.end(),
688 __varcov.begin(), __varcov.end());
689 else
690 _M_init_diagonal(__mean.begin(), __mean.end(),
691 __varcov.begin(), __varcov.end());
692 }
694 std::array<_RealType, _Dimen>
695 mean() const
696 { return _M_mean; }
698 std::array<_RealType, _M_t_size>
699 varcov() const
700 { return _M_t; }
702 friend bool
703 operator==(const param_type& __p1, const param_type& __p2)
704 { return __p1._M_mean == __p2._M_mean && __p1._M_t == __p2._M_t; }
706 private:
707 template <typename _InputIterator1, typename _InputIterator2>
708 void _M_init_full(_InputIterator1 __meanbegin,
709 _InputIterator1 __meanend,
710 _InputIterator2 __varcovbegin,
711 _InputIterator2 __varcovend);
712 template <typename _InputIterator1, typename _InputIterator2>
713 void _M_init_lower(_InputIterator1 __meanbegin,
714 _InputIterator1 __meanend,
715 _InputIterator2 __varcovbegin,
716 _InputIterator2 __varcovend);
717 template <typename _InputIterator1, typename _InputIterator2>
718 void _M_init_diagonal(_InputIterator1 __meanbegin,
719 _InputIterator1 __meanend,
720 _InputIterator2 __varbegin,
721 _InputIterator2 __varend);
723 std::array<_RealType, _Dimen> _M_mean;
724 std::array<_RealType, _M_t_size> _M_t;
725 };
727 public:
728 normal_mv_distribution()
729 : _M_param(), _M_nd()
730 { }
732 template<typename _ForwardIterator1, typename _ForwardIterator2>
733 normal_mv_distribution(_ForwardIterator1 __meanbegin,
734 _ForwardIterator1 __meanend,
735 _ForwardIterator2 __varcovbegin,
736 _ForwardIterator2 __varcovend)
737 : _M_param(__meanbegin, __meanend, __varcovbegin, __varcovend),
738 _M_nd()
739 { }
741 normal_mv_distribution(std::initializer_list<_RealType> __mean,
742 std::initializer_list<_RealType> __varcov)
743 : _M_param(__mean, __varcov), _M_nd()
744 { }
746 explicit
747 normal_mv_distribution(const param_type& __p)
748 : _M_param(__p), _M_nd()
749 { }
751 /**
752 * @brief Resets the distribution state.
753 */
754 void
755 reset()
756 { _M_nd.reset(); }
758 /**
759 * @brief Returns the mean of the distribution.
760 */
761 result_type
762 mean() const
763 { return _M_param.mean(); }
765 /**
766 * @brief Returns the compact form of the variance/covariance
767 * matrix of the distribution.
768 */
769 std::array<_RealType, _Dimen * (_Dimen + 1) / 2>
770 varcov() const
771 { return _M_param.varcov(); }
773 /**
774 * @brief Returns the parameter set of the distribution.
775 */
776 param_type
777 param() const
778 { return _M_param; }
780 /**
781 * @brief Sets the parameter set of the distribution.
782 * @param __param The new parameter set of the distribution.
783 */
784 void
785 param(const param_type& __param)
786 { _M_param = __param; }
788 /**
789 * @brief Returns the greatest lower bound value of the distribution.
790 */
791 result_type
792 min() const
793 { result_type __res;
794 __res.fill(std::numeric_limits<_RealType>::lowest());
795 return __res; }
797 /**
798 * @brief Returns the least upper bound value of the distribution.
799 */
800 result_type
801 max() const
802 { result_type __res;
803 __res.fill(std::numeric_limits<_RealType>::max());
804 return __res; }
806 /**
807 * @brief Generating functions.
808 */
809 template<typename _UniformRandomNumberGenerator>
810 result_type
811 operator()(_UniformRandomNumberGenerator& __urng)
812 { return this->operator()(__urng, _M_param); }
814 template<typename _UniformRandomNumberGenerator>
815 result_type
816 operator()(_UniformRandomNumberGenerator& __urng,
817 const param_type& __p);
819 template<typename _ForwardIterator,
820 typename _UniformRandomNumberGenerator>
821 void
822 __generate(_ForwardIterator __f, _ForwardIterator __t,
823 _UniformRandomNumberGenerator& __urng)
824 { return this->__generate_impl(__f, __t, __urng, _M_param); }
826 template<typename _ForwardIterator,
827 typename _UniformRandomNumberGenerator>
828 void
829 __generate(_ForwardIterator __f, _ForwardIterator __t,
830 _UniformRandomNumberGenerator& __urng,
831 const param_type& __p)
832 { return this->__generate_impl(__f, __t, __urng, __p); }
834 /**
835 * @brief Return true if two multi-variant normal distributions have
836 * the same parameters and the sequences that would
837 * be generated are equal.
838 */
839 template<size_t _Dimen1, typename _RealType1>
840 friend bool
841 operator==(const
842 __gnu_cxx::normal_mv_distribution<_Dimen1, _RealType1>&
843 __d1,
844 const
845 __gnu_cxx::normal_mv_distribution<_Dimen1, _RealType1>&
846 __d2);
848 /**
849 * @brief Inserts a %normal_mv_distribution random number distribution
850 * @p __x into the output stream @p __os.
851 *
852 * @param __os An output stream.
853 * @param __x A %normal_mv_distribution random number distribution.
854 *
855 * @returns The output stream with the state of @p __x inserted or in
856 * an error state.
857 */
858 template<size_t _Dimen1, typename _RealType1,
859 typename _CharT, typename _Traits>
860 friend std::basic_ostream<_CharT, _Traits>&
861 operator<<(std::basic_ostream<_CharT, _Traits>& __os,
862 const
863 __gnu_cxx::normal_mv_distribution<_Dimen1, _RealType1>&
864 __x);
866 /**
867 * @brief Extracts a %normal_mv_distribution random number distribution
868 * @p __x from the input stream @p __is.
869 *
870 * @param __is An input stream.
871 * @param __x A %normal_mv_distribution random number generator engine.
872 *
873 * @returns The input stream with @p __x extracted or in an error
874 * state.
875 */
876 template<size_t _Dimen1, typename _RealType1,
877 typename _CharT, typename _Traits>
878 friend std::basic_istream<_CharT, _Traits>&
879 operator>>(std::basic_istream<_CharT, _Traits>& __is,
880 __gnu_cxx::normal_mv_distribution<_Dimen1, _RealType1>&
881 __x);
883 private:
884 template<typename _ForwardIterator,
885 typename _UniformRandomNumberGenerator>
886 void
887 __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
888 _UniformRandomNumberGenerator& __urng,
889 const param_type& __p);
891 param_type _M_param;
892 std::normal_distribution<_RealType> _M_nd;
893 };
895 /**
896 * @brief Return true if two multi-variate normal distributions are
897 * different.
898 */
899 template<size_t _Dimen, typename _RealType>
900 inline bool
901 operator!=(const __gnu_cxx::normal_mv_distribution<_Dimen, _RealType>&
902 __d1,
903 const __gnu_cxx::normal_mv_distribution<_Dimen, _RealType>&
904 __d2)
905 { return !(__d1 == __d2); }
908 /**
909 * @brief A Rice continuous distribution for random numbers.
910 *
911 * The formula for the Rice probability density function is
912 * @f[
913 * p(x|\nu,\sigma) = \frac{x}{\sigma^2}
914 * \exp\left(-\frac{x^2+\nu^2}{2\sigma^2}\right)
915 * I_0\left(\frac{x \nu}{\sigma^2}\right)
916 * @f]
917 * where @f$I_0(z)@f$ is the modified Bessel function of the first kind
918 * of order 0 and @f$\nu >= 0@f$ and @f$\sigma > 0@f$.
919 *
920 * <table border=1 cellpadding=10 cellspacing=0>
921 * <caption align=top>Distribution Statistics</caption>
922 * <tr><td>Mean</td><td>@f$\sqrt{\pi/2}L_{1/2}(-\nu^2/2\sigma^2)@f$</td></tr>
923 * <tr><td>Variance</td><td>@f$2\sigma^2 + \nu^2
924 * + (\pi\sigma^2/2)L^2_{1/2}(-\nu^2/2\sigma^2)@f$</td></tr>
925 * <tr><td>Range</td><td>@f$[0, \infty)@f$</td></tr>
926 * </table>
927 * where @f$L_{1/2}(x)@f$ is the Laguerre polynomial of order 1/2.
928 */
929 template<typename _RealType = double>
930 class
931 rice_distribution
932 {
933 static_assert(std::is_floating_point<_RealType>::value,
934 "template argument not a floating point type");
935 public:
936 /** The type of the range of the distribution. */
937 typedef _RealType result_type;
938 /** Parameter type. */
939 struct param_type
940 {
941 typedef rice_distribution<result_type> distribution_type;
943 param_type(result_type __nu_val = result_type(0),
944 result_type __sigma_val = result_type(1))
945 : _M_nu(__nu_val), _M_sigma(__sigma_val)
946 {
947 _GLIBCXX_DEBUG_ASSERT(_M_nu >= result_type(0));
948 _GLIBCXX_DEBUG_ASSERT(_M_sigma > result_type(0));
949 }
951 result_type
952 nu() const
953 { return _M_nu; }
955 result_type
956 sigma() const
957 { return _M_sigma; }
959 friend bool
960 operator==(const param_type& __p1, const param_type& __p2)
961 { return __p1._M_nu == __p2._M_nu
962 && __p1._M_sigma == __p2._M_sigma; }
964 private:
965 void _M_initialize();
967 result_type _M_nu;
968 result_type _M_sigma;
969 };
971 /**
972 * @brief Constructors.
973 */
974 explicit
975 rice_distribution(result_type __nu_val = result_type(0),
976 result_type __sigma_val = result_type(1))
977 : _M_param(__nu_val, __sigma_val),
978 _M_ndx(__nu_val, __sigma_val),
979 _M_ndy(result_type(0), __sigma_val)
980 { }
982 explicit
983 rice_distribution(const param_type& __p)
984 : _M_param(__p),
985 _M_ndx(__p.nu(), __p.sigma()),
986 _M_ndy(result_type(0), __p.sigma())
987 { }
989 /**
990 * @brief Resets the distribution state.
991 */
992 void
993 reset()
994 {
995 _M_ndx.reset();
996 _M_ndy.reset();
997 }
999 /**
1000 * @brief Return the parameters of the distribution.
1001 */
1002 result_type
1003 nu() const
1004 { return _M_param.nu(); }
1006 result_type
1007 sigma() const
1008 { return _M_param.sigma(); }
1010 /**
1011 * @brief Returns the parameter set of the distribution.
1012 */
1013 param_type
1014 param() const
1015 { return _M_param; }
1017 /**
1018 * @brief Sets the parameter set of the distribution.
1019 * @param __param The new parameter set of the distribution.
1020 */
1021 void
1022 param(const param_type& __param)
1023 { _M_param = __param; }
1025 /**
1026 * @brief Returns the greatest lower bound value of the distribution.
1027 */
1028 result_type
1029 min() const
1030 { return result_type(0); }
1032 /**
1033 * @brief Returns the least upper bound value of the distribution.
1034 */
1035 result_type
1036 max() const
1037 { return std::numeric_limits<result_type>::max(); }
1039 /**
1040 * @brief Generating functions.
1041 */
1042 template<typename _UniformRandomNumberGenerator>
1043 result_type
1044 operator()(_UniformRandomNumberGenerator& __urng)
1045 {
1046 result_type __x = this->_M_ndx(__urng);
1047 result_type __y = this->_M_ndy(__urng);
1048 #if _GLIBCXX_USE_C99_MATH_TR1
1049 return std::hypot(__x, __y);
1050 #else
1051 return std::sqrt(__x * __x + __y * __y);
1052 #endif
1053 }
1055 template<typename _UniformRandomNumberGenerator>
1056 result_type
1057 operator()(_UniformRandomNumberGenerator& __urng,
1058 const param_type& __p)
1059 {
1060 typename std::normal_distribution<result_type>::param_type
1061 __px(__p.nu(), __p.sigma()), __py(result_type(0), __p.sigma());
1062 result_type __x = this->_M_ndx(__px, __urng);
1063 result_type __y = this->_M_ndy(__py, __urng);
1064 #if _GLIBCXX_USE_C99_MATH_TR1
1065 return std::hypot(__x, __y);
1066 #else
1067 return std::sqrt(__x * __x + __y * __y);
1068 #endif
1069 }
1071 template<typename _ForwardIterator,
1072 typename _UniformRandomNumberGenerator>
1073 void
1074 __generate(_ForwardIterator __f, _ForwardIterator __t,
1075 _UniformRandomNumberGenerator& __urng)
1076 { this->__generate(__f, __t, __urng, _M_param); }
1078 template<typename _ForwardIterator,
1079 typename _UniformRandomNumberGenerator>
1080 void
1081 __generate(_ForwardIterator __f, _ForwardIterator __t,
1082 _UniformRandomNumberGenerator& __urng,
1083 const param_type& __p)
1084 { this->__generate_impl(__f, __t, __urng, __p); }
1086 template<typename _UniformRandomNumberGenerator>
1087 void
1088 __generate(result_type* __f, result_type* __t,
1089 _UniformRandomNumberGenerator& __urng,
1090 const param_type& __p)
1091 { this->__generate_impl(__f, __t, __urng, __p); }
1093 /**
1094 * @brief Return true if two Rice distributions have
1095 * the same parameters and the sequences that would
1096 * be generated are equal.
1097 */
1098 friend bool
1099 operator==(const rice_distribution& __d1,
1100 const rice_distribution& __d2)
1101 { return (__d1._M_param == __d2._M_param
1102 && __d1._M_ndx == __d2._M_ndx
1103 && __d1._M_ndy == __d2._M_ndy); }
1105 /**
1106 * @brief Inserts a %rice_distribution random number distribution
1107 * @p __x into the output stream @p __os.
1108 *
1109 * @param __os An output stream.
1110 * @param __x A %rice_distribution random number distribution.
1111 *
1112 * @returns The output stream with the state of @p __x inserted or in
1113 * an error state.
1114 */
1115 template<typename _RealType1, typename _CharT, typename _Traits>
1116 friend std::basic_ostream<_CharT, _Traits>&
1117 operator<<(std::basic_ostream<_CharT, _Traits>&,
1118 const rice_distribution<_RealType1>&);
1120 /**
1121 * @brief Extracts a %rice_distribution random number distribution
1122 * @p __x from the input stream @p __is.
1123 *
1124 * @param __is An input stream.
1125 * @param __x A %rice_distribution random number
1126 * generator engine.
1127 *
1128 * @returns The input stream with @p __x extracted or in an error state.
1129 */
1130 template<typename _RealType1, typename _CharT, typename _Traits>
1131 friend std::basic_istream<_CharT, _Traits>&
1132 operator>>(std::basic_istream<_CharT, _Traits>&,
1133 rice_distribution<_RealType1>&);
1135 private:
1136 template<typename _ForwardIterator,
1137 typename _UniformRandomNumberGenerator>
1138 void
1139 __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
1140 _UniformRandomNumberGenerator& __urng,
1141 const param_type& __p);
1143 param_type _M_param;
1145 std::normal_distribution<result_type> _M_ndx;
1146 std::normal_distribution<result_type> _M_ndy;
1147 };
1149 /**
1150 * @brief Return true if two Rice distributions are not equal.
1151 */
1152 template<typename _RealType1>
1153 inline bool
1154 operator!=(const rice_distribution<_RealType1>& __d1,
1155 const rice_distribution<_RealType1>& __d2)
1156 { return !(__d1 == __d2); }
1159 /**
1160 * @brief A Nakagami continuous distribution for random numbers.
1161 *
1162 * The formula for the Nakagami probability density function is
1163 * @f[
1164 * p(x|\mu,\omega) = \frac{2\mu^\mu}{\Gamma(\mu)\omega^\mu}
1165 * x^{2\mu-1}e^{-\mu x / \omega}
1166 * @f]
1167 * where @f$\Gamma(z)@f$ is the gamma function and @f$\mu >= 0.5@f$
1168 * and @f$\omega > 0@f$.
1169 */
1170 template<typename _RealType = double>
1171 class
1172 nakagami_distribution
1173 {
1174 static_assert(std::is_floating_point<_RealType>::value,
1175 "template argument not a floating point type");
1177 public:
1178 /** The type of the range of the distribution. */
1179 typedef _RealType result_type;
1180 /** Parameter type. */
1181 struct param_type
1182 {
1183 typedef nakagami_distribution<result_type> distribution_type;
1185 param_type(result_type __mu_val = result_type(1),
1186 result_type __omega_val = result_type(1))
1187 : _M_mu(__mu_val), _M_omega(__omega_val)
1188 {
1189 _GLIBCXX_DEBUG_ASSERT(_M_mu >= result_type(0.5L));
1190 _GLIBCXX_DEBUG_ASSERT(_M_omega > result_type(0));
1191 }
1193 result_type
1194 mu() const
1195 { return _M_mu; }
1197 result_type
1198 omega() const
1199 { return _M_omega; }
1201 friend bool
1202 operator==(const param_type& __p1, const param_type& __p2)
1203 { return __p1._M_mu == __p2._M_mu
1204 && __p1._M_omega == __p2._M_omega; }
1206 private:
1207 void _M_initialize();
1209 result_type _M_mu;
1210 result_type _M_omega;
1211 };
1213 /**
1214 * @brief Constructors.
1215 */
1216 explicit
1217 nakagami_distribution(result_type __mu_val = result_type(1),
1218 result_type __omega_val = result_type(1))
1219 : _M_param(__mu_val, __omega_val),
1220 _M_gd(__mu_val, __omega_val / __mu_val)
1221 { }
1223 explicit
1224 nakagami_distribution(const param_type& __p)
1225 : _M_param(__p),
1226 _M_gd(__p.mu(), __p.omega() / __p.mu())
1227 { }
1229 /**
1230 * @brief Resets the distribution state.
1231 */
1232 void
1233 reset()
1234 { _M_gd.reset(); }
1236 /**
1237 * @brief Return the parameters of the distribution.
1238 */
1239 result_type
1240 mu() const
1241 { return _M_param.mu(); }
1243 result_type
1244 omega() const
1245 { return _M_param.omega(); }
1247 /**
1248 * @brief Returns the parameter set of the distribution.
1249 */
1250 param_type
1251 param() const
1252 { return _M_param; }
1254 /**
1255 * @brief Sets the parameter set of the distribution.
1256 * @param __param The new parameter set of the distribution.
1257 */
1258 void
1259 param(const param_type& __param)
1260 { _M_param = __param; }
1262 /**
1263 * @brief Returns the greatest lower bound value of the distribution.
1264 */
1265 result_type
1266 min() const
1267 { return result_type(0); }
1269 /**
1270 * @brief Returns the least upper bound value of the distribution.
1271 */
1272 result_type
1273 max() const
1274 { return std::numeric_limits<result_type>::max(); }
1276 /**
1277 * @brief Generating functions.
1278 */
1279 template<typename _UniformRandomNumberGenerator>
1280 result_type
1281 operator()(_UniformRandomNumberGenerator& __urng)
1282 { return std::sqrt(this->_M_gd(__urng)); }
1284 template<typename _UniformRandomNumberGenerator>
1285 result_type
1286 operator()(_UniformRandomNumberGenerator& __urng,
1287 const param_type& __p)
1288 {
1289 typename std::gamma_distribution<result_type>::param_type
1290 __pg(__p.mu(), __p.omega() / __p.mu());
1291 return std::sqrt(this->_M_gd(__pg, __urng));
1292 }
1294 template<typename _ForwardIterator,
1295 typename _UniformRandomNumberGenerator>
1296 void
1297 __generate(_ForwardIterator __f, _ForwardIterator __t,
1298 _UniformRandomNumberGenerator& __urng)
1299 { this->__generate(__f, __t, __urng, _M_param); }
1301 template<typename _ForwardIterator,
1302 typename _UniformRandomNumberGenerator>
1303 void
1304 __generate(_ForwardIterator __f, _ForwardIterator __t,
1305 _UniformRandomNumberGenerator& __urng,
1306 const param_type& __p)
1307 { this->__generate_impl(__f, __t, __urng, __p); }
1309 template<typename _UniformRandomNumberGenerator>
1310 void
1311 __generate(result_type* __f, result_type* __t,
1312 _UniformRandomNumberGenerator& __urng,
1313 const param_type& __p)
1314 { this->__generate_impl(__f, __t, __urng, __p); }
1316 /**
1317 * @brief Return true if two Nakagami distributions have
1318 * the same parameters and the sequences that would
1319 * be generated are equal.
1320 */
1321 friend bool
1322 operator==(const nakagami_distribution& __d1,
1323 const nakagami_distribution& __d2)
1324 { return (__d1._M_param == __d2._M_param
1325 && __d1._M_gd == __d2._M_gd); }
1327 /**
1328 * @brief Inserts a %nakagami_distribution random number distribution
1329 * @p __x into the output stream @p __os.
1330 *
1331 * @param __os An output stream.
1332 * @param __x A %nakagami_distribution random number distribution.
1333 *
1334 * @returns The output stream with the state of @p __x inserted or in
1335 * an error state.
1336 */
1337 template<typename _RealType1, typename _CharT, typename _Traits>
1338 friend std::basic_ostream<_CharT, _Traits>&
1339 operator<<(std::basic_ostream<_CharT, _Traits>&,
1340 const nakagami_distribution<_RealType1>&);
1342 /**
1343 * @brief Extracts a %nakagami_distribution random number distribution
1344 * @p __x from the input stream @p __is.
1345 *
1346 * @param __is An input stream.
1347 * @param __x A %nakagami_distribution random number
1348 * generator engine.
1349 *
1350 * @returns The input stream with @p __x extracted or in an error state.
1351 */
1352 template<typename _RealType1, typename _CharT, typename _Traits>
1353 friend std::basic_istream<_CharT, _Traits>&
1354 operator>>(std::basic_istream<_CharT, _Traits>&,
1355 nakagami_distribution<_RealType1>&);
1357 private:
1358 template<typename _ForwardIterator,
1359 typename _UniformRandomNumberGenerator>
1360 void
1361 __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
1362 _UniformRandomNumberGenerator& __urng,
1363 const param_type& __p);
1365 param_type _M_param;
1367 std::gamma_distribution<result_type> _M_gd;
1368 };
1370 /**
1371 * @brief Return true if two Nakagami distributions are not equal.
1372 */
1373 template<typename _RealType>
1374 inline bool
1375 operator!=(const nakagami_distribution<_RealType>& __d1,
1376 const nakagami_distribution<_RealType>& __d2)
1377 { return !(__d1 == __d2); }
1380 /**
1381 * @brief A Pareto continuous distribution for random numbers.
1382 *
1383 * The formula for the Pareto cumulative probability function is
1384 * @f[
1385 * P(x|\alpha,\mu) = 1 - \left(\frac{\mu}{x}\right)^\alpha
1386 * @f]
1387 * The formula for the Pareto probability density function is
1388 * @f[
1389 * p(x|\alpha,\mu) = \frac{\alpha + 1}{\mu}
1390 * \left(\frac{\mu}{x}\right)^{\alpha + 1}
1391 * @f]
1392 * where @f$x >= \mu@f$ and @f$\mu > 0@f$, @f$\alpha > 0@f$.
1393 *
1394 * <table border=1 cellpadding=10 cellspacing=0>
1395 * <caption align=top>Distribution Statistics</caption>
1396 * <tr><td>Mean</td><td>@f$\alpha \mu / (\alpha - 1)@f$
1397 * for @f$\alpha > 1@f$</td></tr>
1398 * <tr><td>Variance</td><td>@f$\alpha \mu^2 / [(\alpha - 1)^2(\alpha - 2)]@f$
1399 * for @f$\alpha > 2@f$</td></tr>
1400 * <tr><td>Range</td><td>@f$[\mu, \infty)@f$</td></tr>
1401 * </table>
1402 */
1403 template<typename _RealType = double>
1404 class
1405 pareto_distribution
1406 {
1407 static_assert(std::is_floating_point<_RealType>::value,
1408 "template argument not a floating point type");
1410 public:
1411 /** The type of the range of the distribution. */
1412 typedef _RealType result_type;
1413 /** Parameter type. */
1414 struct param_type
1415 {
1416 typedef pareto_distribution<result_type> distribution_type;
1418 param_type(result_type __alpha_val = result_type(1),
1419 result_type __mu_val = result_type(1))
1420 : _M_alpha(__alpha_val), _M_mu(__mu_val)
1421 {
1422 _GLIBCXX_DEBUG_ASSERT(_M_alpha > result_type(0));
1423 _GLIBCXX_DEBUG_ASSERT(_M_mu > result_type(0));
1424 }
1426 result_type
1427 alpha() const
1428 { return _M_alpha; }
1430 result_type
1431 mu() const
1432 { return _M_mu; }
1434 friend bool
1435 operator==(const param_type& __p1, const param_type& __p2)
1436 { return __p1._M_alpha == __p2._M_alpha && __p1._M_mu == __p2._M_mu; }
1438 private:
1439 void _M_initialize();
1441 result_type _M_alpha;
1442 result_type _M_mu;
1443 };
1445 /**
1446 * @brief Constructors.
1447 */
1448 explicit
1449 pareto_distribution(result_type __alpha_val = result_type(1),
1450 result_type __mu_val = result_type(1))
1451 : _M_param(__alpha_val, __mu_val),
1452 _M_ud()
1453 { }
1455 explicit
1456 pareto_distribution(const param_type& __p)
1457 : _M_param(__p),
1458 _M_ud()
1459 { }
1461 /**
1462 * @brief Resets the distribution state.
1463 */
1464 void
1465 reset()
1466 {
1467 _M_ud.reset();
1468 }
1470 /**
1471 * @brief Return the parameters of the distribution.
1472 */
1473 result_type
1474 alpha() const
1475 { return _M_param.alpha(); }
1477 result_type
1478 mu() const
1479 { return _M_param.mu(); }
1481 /**
1482 * @brief Returns the parameter set of the distribution.
1483 */
1484 param_type
1485 param() const
1486 { return _M_param; }
1488 /**
1489 * @brief Sets the parameter set of the distribution.
1490 * @param __param The new parameter set of the distribution.
1491 */
1492 void
1493 param(const param_type& __param)
1494 { _M_param = __param; }
1496 /**
1497 * @brief Returns the greatest lower bound value of the distribution.
1498 */
1499 result_type
1500 min() const
1501 { return this->mu(); }
1503 /**
1504 * @brief Returns the least upper bound value of the distribution.
1505 */
1506 result_type
1507 max() const
1508 { return std::numeric_limits<result_type>::max(); }
1510 /**
1511 * @brief Generating functions.
1512 */
1513 template<typename _UniformRandomNumberGenerator>
1514 result_type
1515 operator()(_UniformRandomNumberGenerator& __urng)
1516 {
1517 return this->mu() * std::pow(this->_M_ud(__urng),
1518 -result_type(1) / this->alpha());
1519 }
1521 template<typename _UniformRandomNumberGenerator>
1522 result_type
1523 operator()(_UniformRandomNumberGenerator& __urng,
1524 const param_type& __p)
1525 {
1526 return __p.mu() * std::pow(this->_M_ud(__urng),
1527 -result_type(1) / __p.alpha());
1528 }
1530 template<typename _ForwardIterator,
1531 typename _UniformRandomNumberGenerator>
1532 void
1533 __generate(_ForwardIterator __f, _ForwardIterator __t,
1534 _UniformRandomNumberGenerator& __urng)
1535 { this->__generate(__f, __t, __urng, _M_param); }
1537 template<typename _ForwardIterator,
1538 typename _UniformRandomNumberGenerator>
1539 void
1540 __generate(_ForwardIterator __f, _ForwardIterator __t,
1541 _UniformRandomNumberGenerator& __urng,
1542 const param_type& __p)
1543 { this->__generate_impl(__f, __t, __urng, __p); }
1545 template<typename _UniformRandomNumberGenerator>
1546 void
1547 __generate(result_type* __f, result_type* __t,
1548 _UniformRandomNumberGenerator& __urng,
1549 const param_type& __p)
1550 { this->__generate_impl(__f, __t, __urng, __p); }
1552 /**
1553 * @brief Return true if two Pareto distributions have
1554 * the same parameters and the sequences that would
1555 * be generated are equal.
1556 */
1557 friend bool
1558 operator==(const pareto_distribution& __d1,
1559 const pareto_distribution& __d2)
1560 { return (__d1._M_param == __d2._M_param
1561 && __d1._M_ud == __d2._M_ud); }
1563 /**
1564 * @brief Inserts a %pareto_distribution random number distribution
1565 * @p __x into the output stream @p __os.
1566 *
1567 * @param __os An output stream.
1568 * @param __x A %pareto_distribution random number distribution.
1569 *
1570 * @returns The output stream with the state of @p __x inserted or in
1571 * an error state.
1572 */
1573 template<typename _RealType1, typename _CharT, typename _Traits>
1574 friend std::basic_ostream<_CharT, _Traits>&
1575 operator<<(std::basic_ostream<_CharT, _Traits>&,
1576 const pareto_distribution<_RealType1>&);
1578 /**
1579 * @brief Extracts a %pareto_distribution random number distribution
1580 * @p __x from the input stream @p __is.
1581 *
1582 * @param __is An input stream.
1583 * @param __x A %pareto_distribution random number
1584 * generator engine.
1585 *
1586 * @returns The input stream with @p __x extracted or in an error state.
1587 */
1588 template<typename _RealType1, typename _CharT, typename _Traits>
1589 friend std::basic_istream<_CharT, _Traits>&
1590 operator>>(std::basic_istream<_CharT, _Traits>&,
1591 pareto_distribution<_RealType1>&);
1593 private:
1594 template<typename _ForwardIterator,
1595 typename _UniformRandomNumberGenerator>
1596 void
1597 __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
1598 _UniformRandomNumberGenerator& __urng,
1599 const param_type& __p);
1601 param_type _M_param;
1603 std::uniform_real_distribution<result_type> _M_ud;
1604 };
1606 /**
1607 * @brief Return true if two Pareto distributions are not equal.
1608 */
1609 template<typename _RealType>
1610 inline bool
1611 operator!=(const pareto_distribution<_RealType>& __d1,
1612 const pareto_distribution<_RealType>& __d2)
1613 { return !(__d1 == __d2); }
1616 /**
1617 * @brief A K continuous distribution for random numbers.
1618 *
1619 * The formula for the K probability density function is
1620 * @f[
1621 * p(x|\lambda, \mu, \nu) = \frac{2}{x}
1622 * \left(\frac{\lambda\nu x}{\mu}\right)^{\frac{\lambda + \nu}{2}}
1623 * \frac{1}{\Gamma(\lambda)\Gamma(\nu)}
1624 * K_{\nu - \lambda}\left(2\sqrt{\frac{\lambda\nu x}{\mu}}\right)
1625 * @f]
1626 * where @f$I_0(z)@f$ is the modified Bessel function of the second kind
1627 * of order @f$\nu - \lambda@f$ and @f$\lambda > 0@f$, @f$\mu > 0@f$
1628 * and @f$\nu > 0@f$.
1629 *
1630 * <table border=1 cellpadding=10 cellspacing=0>
1631 * <caption align=top>Distribution Statistics</caption>
1632 * <tr><td>Mean</td><td>@f$\mu@f$</td></tr>
1633 * <tr><td>Variance</td><td>@f$\mu^2\frac{\lambda + \nu + 1}{\lambda\nu}@f$</td></tr>
1634 * <tr><td>Range</td><td>@f$[0, \infty)@f$</td></tr>
1635 * </table>
1636 */
1637 template<typename _RealType = double>
1638 class
1639 k_distribution
1640 {
1641 static_assert(std::is_floating_point<_RealType>::value,
1642 "template argument not a floating point type");
1644 public:
1645 /** The type of the range of the distribution. */
1646 typedef _RealType result_type;
1647 /** Parameter type. */
1648 struct param_type
1649 {
1650 typedef k_distribution<result_type> distribution_type;
1652 param_type(result_type __lambda_val = result_type(1),
1653 result_type __mu_val = result_type(1),
1654 result_type __nu_val = result_type(1))
1655 : _M_lambda(__lambda_val), _M_mu(__mu_val), _M_nu(__nu_val)
1656 {
1657 _GLIBCXX_DEBUG_ASSERT(_M_lambda > result_type(0));
1658 _GLIBCXX_DEBUG_ASSERT(_M_mu > result_type(0));
1659 _GLIBCXX_DEBUG_ASSERT(_M_nu > result_type(0));
1660 }
1662 result_type
1663 lambda() const
1664 { return _M_lambda; }
1666 result_type
1667 mu() const
1668 { return _M_mu; }
1670 result_type
1671 nu() const
1672 { return _M_nu; }
1674 friend bool
1675 operator==(const param_type& __p1, const param_type& __p2)
1676 { return __p1._M_lambda == __p2._M_lambda
1677 && __p1._M_mu == __p2._M_mu
1678 && __p1._M_nu == __p2._M_nu; }
1680 private:
1681 void _M_initialize();
1683 result_type _M_lambda;
1684 result_type _M_mu;
1685 result_type _M_nu;
1686 };
1688 /**
1689 * @brief Constructors.
1690 */
1691 explicit
1692 k_distribution(result_type __lambda_val = result_type(1),
1693 result_type __mu_val = result_type(1),
1694 result_type __nu_val = result_type(1))
1695 : _M_param(__lambda_val, __mu_val, __nu_val),
1696 _M_gd1(__lambda_val, result_type(1) / __lambda_val),
1697 _M_gd2(__nu_val, __mu_val / __nu_val)
1698 { }
1700 explicit
1701 k_distribution(const param_type& __p)
1702 : _M_param(__p),
1703 _M_gd1(__p.lambda(), result_type(1) / __p.lambda()),
1704 _M_gd2(__p.nu(), __p.mu() / __p.nu())
1705 { }
1707 /**
1708 * @brief Resets the distribution state.
1709 */
1710 void
1711 reset()
1712 {
1713 _M_gd1.reset();
1714 _M_gd2.reset();
1715 }
1717 /**
1718 * @brief Return the parameters of the distribution.
1719 */
1720 result_type
1721 lambda() const
1722 { return _M_param.lambda(); }
1724 result_type
1725 mu() const
1726 { return _M_param.mu(); }
1728 result_type
1729 nu() const
1730 { return _M_param.nu(); }
1732 /**
1733 * @brief Returns the parameter set of the distribution.
1734 */
1735 param_type
1736 param() const
1737 { return _M_param; }
1739 /**
1740 * @brief Sets the parameter set of the distribution.
1741 * @param __param The new parameter set of the distribution.
1742 */
1743 void
1744 param(const param_type& __param)
1745 { _M_param = __param; }
1747 /**
1748 * @brief Returns the greatest lower bound value of the distribution.
1749 */
1750 result_type
1751 min() const
1752 { return result_type(0); }
1754 /**
1755 * @brief Returns the least upper bound value of the distribution.
1756 */
1757 result_type
1758 max() const
1759 { return std::numeric_limits<result_type>::max(); }
1761 /**
1762 * @brief Generating functions.
1763 */
1764 template<typename _UniformRandomNumberGenerator>
1765 result_type
1766 operator()(_UniformRandomNumberGenerator&);
1768 template<typename _UniformRandomNumberGenerator>
1769 result_type
1770 operator()(_UniformRandomNumberGenerator&, const param_type&);
1772 template<typename _ForwardIterator,
1773 typename _UniformRandomNumberGenerator>
1774 void
1775 __generate(_ForwardIterator __f, _ForwardIterator __t,
1776 _UniformRandomNumberGenerator& __urng)
1777 { this->__generate(__f, __t, __urng, _M_param); }
1779 template<typename _ForwardIterator,
1780 typename _UniformRandomNumberGenerator>
1781 void
1782 __generate(_ForwardIterator __f, _ForwardIterator __t,
1783 _UniformRandomNumberGenerator& __urng,
1784 const param_type& __p)
1785 { this->__generate_impl(__f, __t, __urng, __p); }
1787 template<typename _UniformRandomNumberGenerator>
1788 void
1789 __generate(result_type* __f, result_type* __t,
1790 _UniformRandomNumberGenerator& __urng,
1791 const param_type& __p)
1792 { this->__generate_impl(__f, __t, __urng, __p); }
1794 /**
1795 * @brief Return true if two K distributions have
1796 * the same parameters and the sequences that would
1797 * be generated are equal.
1798 */
1799 friend bool
1800 operator==(const k_distribution& __d1,
1801 const k_distribution& __d2)
1802 { return (__d1._M_param == __d2._M_param
1803 && __d1._M_gd1 == __d2._M_gd1
1804 && __d1._M_gd2 == __d2._M_gd2); }
1806 /**
1807 * @brief Inserts a %k_distribution random number distribution
1808 * @p __x into the output stream @p __os.
1809 *
1810 * @param __os An output stream.
1811 * @param __x A %k_distribution random number distribution.
1812 *
1813 * @returns The output stream with the state of @p __x inserted or in
1814 * an error state.
1815 */
1816 template<typename _RealType1, typename _CharT, typename _Traits>
1817 friend std::basic_ostream<_CharT, _Traits>&
1818 operator<<(std::basic_ostream<_CharT, _Traits>&,
1819 const k_distribution<_RealType1>&);
1821 /**
1822 * @brief Extracts a %k_distribution random number distribution
1823 * @p __x from the input stream @p __is.
1824 *
1825 * @param __is An input stream.
1826 * @param __x A %k_distribution random number
1827 * generator engine.
1828 *
1829 * @returns The input stream with @p __x extracted or in an error state.
1830 */
1831 template<typename _RealType1, typename _CharT, typename _Traits>
1832 friend std::basic_istream<_CharT, _Traits>&
1833 operator>>(std::basic_istream<_CharT, _Traits>&,
1834 k_distribution<_RealType1>&);
1836 private:
1837 template<typename _ForwardIterator,
1838 typename _UniformRandomNumberGenerator>
1839 void
1840 __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
1841 _UniformRandomNumberGenerator& __urng,
1842 const param_type& __p);
1844 param_type _M_param;
1846 std::gamma_distribution<result_type> _M_gd1;
1847 std::gamma_distribution<result_type> _M_gd2;
1848 };
1850 /**
1851 * @brief Return true if two K distributions are not equal.
1852 */
1853 template<typename _RealType>
1854 inline bool
1855 operator!=(const k_distribution<_RealType>& __d1,
1856 const k_distribution<_RealType>& __d2)
1857 { return !(__d1 == __d2); }
1860 /**
1861 * @brief An arcsine continuous distribution for random numbers.
1862 *
1863 * The formula for the arcsine probability density function is
1864 * @f[
1865 * p(x|a,b) = \frac{1}{\pi \sqrt{(x - a)(b - x)}}
1866 * @f]
1867 * where @f$x >= a@f$ and @f$x <= b@f$.
1868 *
1869 * <table border=1 cellpadding=10 cellspacing=0>
1870 * <caption align=top>Distribution Statistics</caption>
1871 * <tr><td>Mean</td><td>@f$ (a + b) / 2 @f$</td></tr>
1872 * <tr><td>Variance</td><td>@f$ (b - a)^2 / 8 @f$</td></tr>
1873 * <tr><td>Range</td><td>@f$[a, b]@f$</td></tr>
1874 * </table>
1875 */
1876 template<typename _RealType = double>
1877 class
1878 arcsine_distribution
1879 {
1880 static_assert(std::is_floating_point<_RealType>::value,
1881 "template argument not a floating point type");
1883 public:
1884 /** The type of the range of the distribution. */
1885 typedef _RealType result_type;
1886 /** Parameter type. */
1887 struct param_type
1888 {
1889 typedef arcsine_distribution<result_type> distribution_type;
1891 param_type(result_type __a = result_type(0),
1892 result_type __b = result_type(1))
1893 : _M_a(__a), _M_b(__b)
1894 {
1895 _GLIBCXX_DEBUG_ASSERT(_M_a <= _M_b);
1896 }
1898 result_type
1899 a() const
1900 { return _M_a; }
1902 result_type
1903 b() const
1904 { return _M_b; }
1906 friend bool
1907 operator==(const param_type& __p1, const param_type& __p2)
1908 { return __p1._M_a == __p2._M_a && __p1._M_b == __p2._M_b; }
1910 private:
1911 void _M_initialize();
1913 result_type _M_a;
1914 result_type _M_b;
1915 };
1917 /**
1918 * @brief Constructors.
1919 */
1920 explicit
1921 arcsine_distribution(result_type __a = result_type(0),
1922 result_type __b = result_type(1))
1923 : _M_param(__a, __b),
1924 _M_ud(-1.5707963267948966192313216916397514L,
1925 +1.5707963267948966192313216916397514L)
1926 { }
1928 explicit
1929 arcsine_distribution(const param_type& __p)
1930 : _M_param(__p),
1931 _M_ud(-1.5707963267948966192313216916397514L,
1932 +1.5707963267948966192313216916397514L)
1933 { }
1935 /**
1936 * @brief Resets the distribution state.
1937 */
1938 void
1939 reset()
1940 { _M_ud.reset(); }
1942 /**
1943 * @brief Return the parameters of the distribution.
1944 */
1945 result_type
1946 a() const
1947 { return _M_param.a(); }
1949 result_type
1950 b() const
1951 { return _M_param.b(); }
1953 /**
1954 * @brief Returns the parameter set of the distribution.
1955 */
1956 param_type
1957 param() const
1958 { return _M_param; }
1960 /**
1961 * @brief Sets the parameter set of the distribution.
1962 * @param __param The new parameter set of the distribution.
1963 */
1964 void
1965 param(const param_type& __param)
1966 { _M_param = __param; }
1968 /**
1969 * @brief Returns the greatest lower bound value of the distribution.
1970 */
1971 result_type
1972 min() const
1973 { return this->a(); }
1975 /**
1976 * @brief Returns the least upper bound value of the distribution.
1977 */
1978 result_type
1979 max() const
1980 { return this->b(); }
1982 /**
1983 * @brief Generating functions.
1984 */
1985 template<typename _UniformRandomNumberGenerator>
1986 result_type
1987 operator()(_UniformRandomNumberGenerator& __urng)
1988 {
1989 result_type __x = std::sin(this->_M_ud(__urng));
1990 return (__x * (this->b() - this->a())
1991 + this->a() + this->b()) / result_type(2);
1992 }
1994 template<typename _UniformRandomNumberGenerator>
1995 result_type
1996 operator()(_UniformRandomNumberGenerator& __urng,
1997 const param_type& __p)
1998 {
1999 result_type __x = std::sin(this->_M_ud(__urng));
2000 return (__x * (__p.b() - __p.a())
2001 + __p.a() + __p.b()) / result_type(2);
2002 }
2004 template<typename _ForwardIterator,
2005 typename _UniformRandomNumberGenerator>
2006 void
2007 __generate(_ForwardIterator __f, _ForwardIterator __t,
2008 _UniformRandomNumberGenerator& __urng)
2009 { this->__generate(__f, __t, __urng, _M_param); }
2011 template<typename _ForwardIterator,
2012 typename _UniformRandomNumberGenerator>
2013 void
2014 __generate(_ForwardIterator __f, _ForwardIterator __t,
2015 _UniformRandomNumberGenerator& __urng,
2016 const param_type& __p)
2017 { this->__generate_impl(__f, __t, __urng, __p); }
2019 template<typename _UniformRandomNumberGenerator>
2020 void
2021 __generate(result_type* __f, result_type* __t,
2022 _UniformRandomNumberGenerator& __urng,
2023 const param_type& __p)
2024 { this->__generate_impl(__f, __t, __urng, __p); }
2026 /**
2027 * @brief Return true if two arcsine distributions have
2028 * the same parameters and the sequences that would
2029 * be generated are equal.
2030 */
2031 friend bool
2032 operator==(const arcsine_distribution& __d1,
2033 const arcsine_distribution& __d2)
2034 { return (__d1._M_param == __d2._M_param
2035 && __d1._M_ud == __d2._M_ud); }
2037 /**
2038 * @brief Inserts a %arcsine_distribution random number distribution
2039 * @p __x into the output stream @p __os.
2040 *
2041 * @param __os An output stream.
2042 * @param __x A %arcsine_distribution random number distribution.
2043 *
2044 * @returns The output stream with the state of @p __x inserted or in
2045 * an error state.
2046 */
2047 template<typename _RealType1, typename _CharT, typename _Traits>
2048 friend std::basic_ostream<_CharT, _Traits>&
2049 operator<<(std::basic_ostream<_CharT, _Traits>&,
2050 const arcsine_distribution<_RealType1>&);
2052 /**
2053 * @brief Extracts a %arcsine_distribution random number distribution
2054 * @p __x from the input stream @p __is.
2055 *
2056 * @param __is An input stream.
2057 * @param __x A %arcsine_distribution random number
2058 * generator engine.
2059 *
2060 * @returns The input stream with @p __x extracted or in an error state.
2061 */
2062 template<typename _RealType1, typename _CharT, typename _Traits>
2063 friend std::basic_istream<_CharT, _Traits>&
2064 operator>>(std::basic_istream<_CharT, _Traits>&,
2065 arcsine_distribution<_RealType1>&);
2067 private:
2068 template<typename _ForwardIterator,
2069 typename _UniformRandomNumberGenerator>
2070 void
2071 __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
2072 _UniformRandomNumberGenerator& __urng,
2073 const param_type& __p);
2075 param_type _M_param;
2077 std::uniform_real_distribution<result_type> _M_ud;
2078 };
2080 /**
2081 * @brief Return true if two arcsine distributions are not equal.
2082 */
2083 template<typename _RealType>
2084 inline bool
2085 operator!=(const arcsine_distribution<_RealType>& __d1,
2086 const arcsine_distribution<_RealType>& __d2)
2087 { return !(__d1 == __d2); }
2090 /**
2091 * @brief A Hoyt continuous distribution for random numbers.
2092 *
2093 * The formula for the Hoyt probability density function is
2094 * @f[
2095 * p(x|q,\omega) = \frac{(1 + q^2)x}{q\omega}
2096 * \exp\left(-\frac{(1 + q^2)^2 x^2}{4 q^2 \omega}\right)
2097 * I_0\left(\frac{(1 - q^4) x^2}{4 q^2 \omega}\right)
2098 * @f]
2099 * where @f$I_0(z)@f$ is the modified Bessel function of the first kind
2100 * of order 0 and @f$0 < q < 1@f$.
2101 *
2102 * <table border=1 cellpadding=10 cellspacing=0>
2103 * <caption align=top>Distribution Statistics</caption>
2104 * <tr><td>Mean</td><td>@f$ \sqrt{\frac{2}{\pi}} \sqrt{\frac{\omega}{1 + q^2}}
2105 * E(1 - q^2) @f$</td></tr>
2106 * <tr><td>Variance</td><td>@f$ \omega \left(1 - \frac{2E^2(1 - q^2)}
2107 * {\pi (1 + q^2)}\right) @f$</td></tr>
2108 * <tr><td>Range</td><td>@f$[0, \infty)@f$</td></tr>
2109 * </table>
2110 * where @f$E(x)@f$ is the elliptic function of the second kind.
2111 */
2112 template<typename _RealType = double>
2113 class
2114 hoyt_distribution
2115 {
2116 static_assert(std::is_floating_point<_RealType>::value,
2117 "template argument not a floating point type");
2119 public:
2120 /** The type of the range of the distribution. */
2121 typedef _RealType result_type;
2122 /** Parameter type. */
2123 struct param_type
2124 {
2125 typedef hoyt_distribution<result_type> distribution_type;
2127 param_type(result_type __q = result_type(0.5L),
2128 result_type __omega = result_type(1))
2129 : _M_q(__q), _M_omega(__omega)
2130 {
2131 _GLIBCXX_DEBUG_ASSERT(_M_q > result_type(0));
2132 _GLIBCXX_DEBUG_ASSERT(_M_q < result_type(1));
2133 }
2135 result_type
2136 q() const
2137 { return _M_q; }
2139 result_type
2140 omega() const
2141 { return _M_omega; }
2143 friend bool
2144 operator==(const param_type& __p1, const param_type& __p2)
2145 { return __p1._M_q == __p2._M_q
2146 && __p1._M_omega == __p2._M_omega; }
2148 private:
2149 void _M_initialize();
2151 result_type _M_q;
2152 result_type _M_omega;
2153 };
2155 /**
2156 * @brief Constructors.
2157 */
2158 explicit
2159 hoyt_distribution(result_type __q = result_type(0.5L),
2160 result_type __omega = result_type(1))
2161 : _M_param(__q, __omega),
2162 _M_ad(result_type(0.5L) * (result_type(1) + __q * __q),
2163 result_type(0.5L) * (result_type(1) + __q * __q)
2164 / (__q * __q)),
2165 _M_ed(result_type(1))
2166 { }
2168 explicit
2169 hoyt_distribution(const param_type& __p)
2170 : _M_param(__p),
2171 _M_ad(result_type(0.5L) * (result_type(1) + __p.q() * __p.q()),
2172 result_type(0.5L) * (result_type(1) + __p.q() * __p.q())
2173 / (__p.q() * __p.q())),
2174 _M_ed(result_type(1))
2175 { }
2177 /**
2178 * @brief Resets the distribution state.
2179 */
2180 void
2181 reset()
2182 {
2183 _M_ad.reset();
2184 _M_ed.reset();
2185 }
2187 /**
2188 * @brief Return the parameters of the distribution.
2189 */
2190 result_type
2191 q() const
2192 { return _M_param.q(); }
2194 result_type
2195 omega() const
2196 { return _M_param.omega(); }
2198 /**
2199 * @brief Returns the parameter set of the distribution.
2200 */
2201 param_type
2202 param() const
2203 { return _M_param; }
2205 /**
2206 * @brief Sets the parameter set of the distribution.
2207 * @param __param The new parameter set of the distribution.
2208 */
2209 void
2210 param(const param_type& __param)
2211 { _M_param = __param; }
2213 /**
2214 * @brief Returns the greatest lower bound value of the distribution.
2215 */
2216 result_type
2217 min() const
2218 { return result_type(0); }
2220 /**
2221 * @brief Returns the least upper bound value of the distribution.
2222 */
2223 result_type
2224 max() const
2225 { return std::numeric_limits<result_type>::max(); }
2227 /**
2228 * @brief Generating functions.
2229 */
2230 template<typename _UniformRandomNumberGenerator>
2231 result_type
2232 operator()(_UniformRandomNumberGenerator& __urng);
2234 template<typename _UniformRandomNumberGenerator>
2235 result_type
2236 operator()(_UniformRandomNumberGenerator& __urng,
2237 const param_type& __p);
2239 template<typename _ForwardIterator,
2240 typename _UniformRandomNumberGenerator>
2241 void
2242 __generate(_ForwardIterator __f, _ForwardIterator __t,
2243 _UniformRandomNumberGenerator& __urng)
2244 { this->__generate(__f, __t, __urng, _M_param); }
2246 template<typename _ForwardIterator,
2247 typename _UniformRandomNumberGenerator>
2248 void
2249 __generate(_ForwardIterator __f, _ForwardIterator __t,
2250 _UniformRandomNumberGenerator& __urng,
2251 const param_type& __p)
2252 { this->__generate_impl(__f, __t, __urng, __p); }
2254 template<typename _UniformRandomNumberGenerator>
2255 void
2256 __generate(result_type* __f, result_type* __t,
2257 _UniformRandomNumberGenerator& __urng,
2258 const param_type& __p)
2259 { this->__generate_impl(__f, __t, __urng, __p); }
2261 /**
2262 * @brief Return true if two Hoyt distributions have
2263 * the same parameters and the sequences that would
2264 * be generated are equal.
2265 */
2266 friend bool
2267 operator==(const hoyt_distribution& __d1,
2268 const hoyt_distribution& __d2)
2269 { return (__d1._M_param == __d2._M_param
2270 && __d1._M_ad == __d2._M_ad
2271 && __d1._M_ed == __d2._M_ed); }
2273 /**
2274 * @brief Inserts a %hoyt_distribution random number distribution
2275 * @p __x into the output stream @p __os.
2276 *
2277 * @param __os An output stream.
2278 * @param __x A %hoyt_distribution random number distribution.
2279 *
2280 * @returns The output stream with the state of @p __x inserted or in
2281 * an error state.
2282 */
2283 template<typename _RealType1, typename _CharT, typename _Traits>
2284 friend std::basic_ostream<_CharT, _Traits>&
2285 operator<<(std::basic_ostream<_CharT, _Traits>&,
2286 const hoyt_distribution<_RealType1>&);
2288 /**
2289 * @brief Extracts a %hoyt_distribution random number distribution
2290 * @p __x from the input stream @p __is.
2291 *
2292 * @param __is An input stream.
2293 * @param __x A %hoyt_distribution random number
2294 * generator engine.
2295 *
2296 * @returns The input stream with @p __x extracted or in an error state.
2297 */
2298 template<typename _RealType1, typename _CharT, typename _Traits>
2299 friend std::basic_istream<_CharT, _Traits>&
2300 operator>>(std::basic_istream<_CharT, _Traits>&,
2301 hoyt_distribution<_RealType1>&);
2303 private:
2304 template<typename _ForwardIterator,
2305 typename _UniformRandomNumberGenerator>
2306 void
2307 __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
2308 _UniformRandomNumberGenerator& __urng,
2309 const param_type& __p);
2311 param_type _M_param;
2313 __gnu_cxx::arcsine_distribution<result_type> _M_ad;
2314 std::exponential_distribution<result_type> _M_ed;
2315 };
2317 /**
2318 * @brief Return true if two Hoyt distributions are not equal.
2319 */
2320 template<typename _RealType>
2321 inline bool
2322 operator!=(const hoyt_distribution<_RealType>& __d1,
2323 const hoyt_distribution<_RealType>& __d2)
2324 { return !(__d1 == __d2); }
2327 /**
2328 * @brief A triangular distribution for random numbers.
2329 *
2330 * The formula for the triangular probability density function is
2331 * @f[
2332 * / 0 for x < a
2333 * p(x|a,b,c) = | \frac{2(x-a)}{(c-a)(b-a)} for a <= x <= b
2334 * | \frac{2(c-x)}{(c-a)(c-b)} for b < x <= c
2335 * \ 0 for c < x
2336 * @f]
2337 *
2338 * <table border=1 cellpadding=10 cellspacing=0>
2339 * <caption align=top>Distribution Statistics</caption>
2340 * <tr><td>Mean</td><td>@f$ \frac{a+b+c}{2} @f$</td></tr>
2341 * <tr><td>Variance</td><td>@f$ \frac{a^2+b^2+c^2-ab-ac-bc}
2342 * {18}@f$</td></tr>
2343 * <tr><td>Range</td><td>@f$[a, c]@f$</td></tr>
2344 * </table>
2345 */
2346 template<typename _RealType = double>
2347 class triangular_distribution
2348 {
2349 static_assert(std::is_floating_point<_RealType>::value,
2350 "template argument not a floating point type");
2352 public:
2353 /** The type of the range of the distribution. */
2354 typedef _RealType result_type;
2355 /** Parameter type. */
2356 struct param_type
2357 {
2358 friend class triangular_distribution<_RealType>;
2360 explicit
2361 param_type(_RealType __a = _RealType(0),
2362 _RealType __b = _RealType(0.5),
2363 _RealType __c = _RealType(1))
2364 : _M_a(__a), _M_b(__b), _M_c(__c)
2365 {
2366 _GLIBCXX_DEBUG_ASSERT(_M_a <= _M_b);
2367 _GLIBCXX_DEBUG_ASSERT(_M_b <= _M_c);
2368 _GLIBCXX_DEBUG_ASSERT(_M_a < _M_c);
2370 _M_r_ab = (_M_b - _M_a) / (_M_c - _M_a);
2371 _M_f_ab_ac = (_M_b - _M_a) * (_M_c - _M_a);
2372 _M_f_bc_ac = (_M_c - _M_b) * (_M_c - _M_a);
2373 }
2375 _RealType
2376 a() const
2377 { return _M_a; }
2379 _RealType
2380 b() const
2381 { return _M_b; }
2383 _RealType
2384 c() const
2385 { return _M_c; }
2387 friend bool
2388 operator==(const param_type& __p1, const param_type& __p2)
2389 { return (__p1._M_a == __p2._M_a && __p1._M_b == __p2._M_b
2390 && __p1._M_c == __p2._M_c); }
2392 private:
2394 _RealType _M_a;
2395 _RealType _M_b;
2396 _RealType _M_c;
2397 _RealType _M_r_ab;
2398 _RealType _M_f_ab_ac;
2399 _RealType _M_f_bc_ac;
2400 };
2402 /**
2403 * @brief Constructs a triangle distribution with parameters
2404 * @f$ a @f$, @f$ b @f$ and @f$ c @f$.
2405 */
2406 explicit
2407 triangular_distribution(result_type __a = result_type(0),
2408 result_type __b = result_type(0.5),
2409 result_type __c = result_type(1))
2410 : _M_param(__a, __b, __c)
2411 { }
2413 explicit
2414 triangular_distribution(const param_type& __p)
2415 : _M_param(__p)
2416 { }
2418 /**
2419 * @brief Resets the distribution state.
2420 */
2421 void
2422 reset()
2423 { }
2425 /**
2426 * @brief Returns the @f$ a @f$ of the distribution.
2427 */
2428 result_type
2429 a() const
2430 { return _M_param.a(); }
2432 /**
2433 * @brief Returns the @f$ b @f$ of the distribution.
2434 */
2435 result_type
2436 b() const
2437 { return _M_param.b(); }
2439 /**
2440 * @brief Returns the @f$ c @f$ of the distribution.
2441 */
2442 result_type
2443 c() const
2444 { return _M_param.c(); }
2446 /**
2447 * @brief Returns the parameter set of the distribution.
2448 */
2449 param_type
2450 param() const
2451 { return _M_param; }
2453 /**
2454 * @brief Sets the parameter set of the distribution.
2455 * @param __param The new parameter set of the distribution.
2456 */
2457 void
2458 param(const param_type& __param)
2459 { _M_param = __param; }
2461 /**
2462 * @brief Returns the greatest lower bound value of the distribution.
2463 */
2464 result_type
2465 min() const
2466 { return _M_param._M_a; }
2468 /**
2469 * @brief Returns the least upper bound value of the distribution.
2470 */
2471 result_type
2472 max() const
2473 { return _M_param._M_c; }
2475 /**
2476 * @brief Generating functions.
2477 */
2478 template<typename _UniformRandomNumberGenerator>
2479 result_type
2480 operator()(_UniformRandomNumberGenerator& __urng)
2481 { return this->operator()(__urng, _M_param); }
2483 template<typename _UniformRandomNumberGenerator>
2484 result_type
2485 operator()(_UniformRandomNumberGenerator& __urng,
2486 const param_type& __p)
2487 {
2488 std::__detail::_Adaptor<_UniformRandomNumberGenerator, result_type>
2489 __aurng(__urng);
2490 result_type __rnd = __aurng();
2491 if (__rnd <= __p._M_r_ab)
2492 return __p.a() + std::sqrt(__rnd * __p._M_f_ab_ac);
2493 else
2494 return __p.c() - std::sqrt((result_type(1) - __rnd)
2495 * __p._M_f_bc_ac);
2496 }
2498 template<typename _ForwardIterator,
2499 typename _UniformRandomNumberGenerator>
2500 void
2501 __generate(_ForwardIterator __f, _ForwardIterator __t,
2502 _UniformRandomNumberGenerator& __urng)
2503 { this->__generate(__f, __t, __urng, _M_param); }
2505 template<typename _ForwardIterator,
2506 typename _UniformRandomNumberGenerator>
2507 void
2508 __generate(_ForwardIterator __f, _ForwardIterator __t,
2509 _UniformRandomNumberGenerator& __urng,
2510 const param_type& __p)
2511 { this->__generate_impl(__f, __t, __urng, __p); }
2513 template<typename _UniformRandomNumberGenerator>
2514 void
2515 __generate(result_type* __f, result_type* __t,
2516 _UniformRandomNumberGenerator& __urng,
2517 const param_type& __p)
2518 { this->__generate_impl(__f, __t, __urng, __p); }
2520 /**
2521 * @brief Return true if two triangle distributions have the same
2522 * parameters and the sequences that would be generated
2523 * are equal.
2524 */
2525 friend bool
2526 operator==(const triangular_distribution& __d1,
2527 const triangular_distribution& __d2)
2528 { return __d1._M_param == __d2._M_param; }
2530 /**
2531 * @brief Inserts a %triangular_distribution random number distribution
2532 * @p __x into the output stream @p __os.
2533 *
2534 * @param __os An output stream.
2535 * @param __x A %triangular_distribution random number distribution.
2536 *
2537 * @returns The output stream with the state of @p __x inserted or in
2538 * an error state.
2539 */
2540 template<typename _RealType1, typename _CharT, typename _Traits>
2541 friend std::basic_ostream<_CharT, _Traits>&
2542 operator<<(std::basic_ostream<_CharT, _Traits>& __os,
2543 const __gnu_cxx::triangular_distribution<_RealType1>& __x);
2545 /**
2546 * @brief Extracts a %triangular_distribution random number distribution
2547 * @p __x from the input stream @p __is.
2548 *
2549 * @param __is An input stream.
2550 * @param __x A %triangular_distribution random number generator engine.
2551 *
2552 * @returns The input stream with @p __x extracted or in an error state.
2553 */
2554 template<typename _RealType1, typename _CharT, typename _Traits>
2555 friend std::basic_istream<_CharT, _Traits>&
2556 operator>>(std::basic_istream<_CharT, _Traits>& __is,
2557 __gnu_cxx::triangular_distribution<_RealType1>& __x);
2559 private:
2560 template<typename _ForwardIterator,
2561 typename _UniformRandomNumberGenerator>
2562 void
2563 __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
2564 _UniformRandomNumberGenerator& __urng,
2565 const param_type& __p);
2567 param_type _M_param;
2568 };
2570 /**
2571 * @brief Return true if two triangle distributions are different.
2572 */
2573 template<typename _RealType>
2574 inline bool
2575 operator!=(const __gnu_cxx::triangular_distribution<_RealType>& __d1,
2576 const __gnu_cxx::triangular_distribution<_RealType>& __d2)
2577 { return !(__d1 == __d2); }
2580 /**
2581 * @brief A von Mises distribution for random numbers.
2582 *
2583 * The formula for the von Mises probability density function is
2584 * @f[
2585 * p(x|\mu,\kappa) = \frac{e^{\kappa \cos(x-\mu)}}
2586 * {2\pi I_0(\kappa)}
2587 * @f]
2588 *
2589 * The generating functions use the method according to:
2590 *
2591 * D. J. Best and N. I. Fisher, 1979. "Efficient Simulation of the
2592 * von Mises Distribution", Journal of the Royal Statistical Society.
2593 * Series C (Applied Statistics), Vol. 28, No. 2, pp. 152-157.
2594 *
2595 * <table border=1 cellpadding=10 cellspacing=0>
2596 * <caption align=top>Distribution Statistics</caption>
2597 * <tr><td>Mean</td><td>@f$ \mu @f$</td></tr>
2598 * <tr><td>Variance</td><td>@f$ 1-I_1(\kappa)/I_0(\kappa) @f$</td></tr>
2599 * <tr><td>Range</td><td>@f$[-\pi, \pi]@f$</td></tr>
2600 * </table>
2601 */
2602 template<typename _RealType = double>
2603 class von_mises_distribution
2604 {
2605 static_assert(std::is_floating_point<_RealType>::value,
2606 "template argument not a floating point type");
2608 public:
2609 /** The type of the range of the distribution. */
2610 typedef _RealType result_type;
2611 /** Parameter type. */
2612 struct param_type
2613 {
2614 friend class von_mises_distribution<_RealType>;
2616 explicit
2617 param_type(_RealType __mu = _RealType(0),
2618 _RealType __kappa = _RealType(1))
2619 : _M_mu(__mu), _M_kappa(__kappa)
2620 {
2621 const _RealType __pi = __gnu_cxx::__math_constants<_RealType>::__pi;
2622 _GLIBCXX_DEBUG_ASSERT(_M_mu >= -__pi && _M_mu <= __pi);
2623 _GLIBCXX_DEBUG_ASSERT(_M_kappa >= _RealType(0));
2625 auto __tau = std::sqrt(_RealType(4) * _M_kappa * _M_kappa
2626 + _RealType(1)) + _RealType(1);
2627 auto __rho = ((__tau - std::sqrt(_RealType(2) * __tau))
2628 / (_RealType(2) * _M_kappa));
2629 _M_r = (_RealType(1) + __rho * __rho) / (_RealType(2) * __rho);
2630 }
2632 _RealType
2633 mu() const
2634 { return _M_mu; }
2636 _RealType
2637 kappa() const
2638 { return _M_kappa; }
2640 friend bool
2641 operator==(const param_type& __p1, const param_type& __p2)
2642 { return (__p1._M_mu == __p2._M_mu
2643 && __p1._M_kappa == __p2._M_kappa); }
2645 private:
2646 _RealType _M_mu;
2647 _RealType _M_kappa;
2648 _RealType _M_r;
2649 };
2651 /**
2652 * @brief Constructs a von Mises distribution with parameters
2653 * @f$\mu@f$ and @f$\kappa@f$.
2654 */
2655 explicit
2656 von_mises_distribution(result_type __mu = result_type(0),
2657 result_type __kappa = result_type(1))
2658 : _M_param(__mu, __kappa)
2659 { }
2661 explicit
2662 von_mises_distribution(const param_type& __p)
2663 : _M_param(__p)
2664 { }
2666 /**
2667 * @brief Resets the distribution state.
2668 */
2669 void
2670 reset()
2671 { }
2673 /**
2674 * @brief Returns the @f$ \mu @f$ of the distribution.
2675 */
2676 result_type
2677 mu() const
2678 { return _M_param.mu(); }
2680 /**
2681 * @brief Returns the @f$ \kappa @f$ of the distribution.
2682 */
2683 result_type
2684 kappa() const
2685 { return _M_param.kappa(); }
2687 /**
2688 * @brief Returns the parameter set of the distribution.
2689 */
2690 param_type
2691 param() const
2692 { return _M_param; }
2694 /**
2695 * @brief Sets the parameter set of the distribution.
2696 * @param __param The new parameter set of the distribution.
2697 */
2698 void
2699 param(const param_type& __param)
2700 { _M_param = __param; }
2702 /**
2703 * @brief Returns the greatest lower bound value of the distribution.
2704 */
2705 result_type
2706 min() const
2707 {
2708 return -__gnu_cxx::__math_constants<result_type>::__pi;
2709 }
2711 /**
2712 * @brief Returns the least upper bound value of the distribution.
2713 */
2714 result_type
2715 max() const
2716 {
2717 return __gnu_cxx::__math_constants<result_type>::__pi;
2718 }
2720 /**
2721 * @brief Generating functions.
2722 */
2723 template<typename _UniformRandomNumberGenerator>
2724 result_type
2725 operator()(_UniformRandomNumberGenerator& __urng)
2726 { return this->operator()(__urng, _M_param); }
2728 template<typename _UniformRandomNumberGenerator>
2729 result_type
2730 operator()(_UniformRandomNumberGenerator& __urng,
2731 const param_type& __p);
2733 template<typename _ForwardIterator,
2734 typename _UniformRandomNumberGenerator>
2735 void
2736 __generate(_ForwardIterator __f, _ForwardIterator __t,
2737 _UniformRandomNumberGenerator& __urng)
2738 { this->__generate(__f, __t, __urng, _M_param); }
2740 template<typename _ForwardIterator,
2741 typename _UniformRandomNumberGenerator>
2742 void
2743 __generate(_ForwardIterator __f, _ForwardIterator __t,
2744 _UniformRandomNumberGenerator& __urng,
2745 const param_type& __p)
2746 { this->__generate_impl(__f, __t, __urng, __p); }
2748 template<typename _UniformRandomNumberGenerator>
2749 void
2750 __generate(result_type* __f, result_type* __t,
2751 _UniformRandomNumberGenerator& __urng,
2752 const param_type& __p)
2753 { this->__generate_impl(__f, __t, __urng, __p); }
2755 /**
2756 * @brief Return true if two von Mises distributions have the same
2757 * parameters and the sequences that would be generated
2758 * are equal.
2759 */
2760 friend bool
2761 operator==(const von_mises_distribution& __d1,
2762 const von_mises_distribution& __d2)
2763 { return __d1._M_param == __d2._M_param; }
2765 /**
2766 * @brief Inserts a %von_mises_distribution random number distribution
2767 * @p __x into the output stream @p __os.
2768 *
2769 * @param __os An output stream.
2770 * @param __x A %von_mises_distribution random number distribution.
2771 *
2772 * @returns The output stream with the state of @p __x inserted or in
2773 * an error state.
2774 */
2775 template<typename _RealType1, typename _CharT, typename _Traits>
2776 friend std::basic_ostream<_CharT, _Traits>&
2777 operator<<(std::basic_ostream<_CharT, _Traits>& __os,
2778 const __gnu_cxx::von_mises_distribution<_RealType1>& __x);
2780 /**
2781 * @brief Extracts a %von_mises_distribution random number distribution
2782 * @p __x from the input stream @p __is.
2783 *
2784 * @param __is An input stream.
2785 * @param __x A %von_mises_distribution random number generator engine.
2786 *
2787 * @returns The input stream with @p __x extracted or in an error state.
2788 */
2789 template<typename _RealType1, typename _CharT, typename _Traits>
2790 friend std::basic_istream<_CharT, _Traits>&
2791 operator>>(std::basic_istream<_CharT, _Traits>& __is,
2792 __gnu_cxx::von_mises_distribution<_RealType1>& __x);
2794 private:
2795 template<typename _ForwardIterator,
2796 typename _UniformRandomNumberGenerator>
2797 void
2798 __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
2799 _UniformRandomNumberGenerator& __urng,
2800 const param_type& __p);
2802 param_type _M_param;
2803 };
2805 /**
2806 * @brief Return true if two von Mises distributions are different.
2807 */
2808 template<typename _RealType>
2809 inline bool
2810 operator!=(const __gnu_cxx::von_mises_distribution<_RealType>& __d1,
2811 const __gnu_cxx::von_mises_distribution<_RealType>& __d2)
2812 { return !(__d1 == __d2); }
2815 /**
2816 * @brief A discrete hypergeometric random number distribution.
2817 *
2818 * The hypergeometric distribution is a discrete probability distribution
2819 * that describes the probability of @p k successes in @p n draws @a without
2820 * replacement from a finite population of size @p N containing exactly @p K
2821 * successes.
2822 *
2823 * The formula for the hypergeometric probability density function is
2824 * @f[
2825 * p(k|N,K,n) = \frac{\binom{K}{k} \binom{N-K}{n-k}}{\binom{N}{n}}
2826 * @f]
2827 * where @f$N@f$ is the total population of the distribution,
2828 * @f$K@f$ is the total population of the distribution.
2829 *
2830 * <table border=1 cellpadding=10 cellspacing=0>
2831 * <caption align=top>Distribution Statistics</caption>
2832 * <tr><td>Mean</td><td>@f$ n\frac{K}{N} @f$</td></tr>
2833 * <tr><td>Variance</td><td>@f$ n\frac{K}{N}\frac{N-K}{N}\frac{N-n}{N-1}
2834 * @f$</td></tr>
2835 * <tr><td>Range</td><td>@f$[max(0, n+K-N), min(K, n)]@f$</td></tr>
2836 * </table>
2837 */
2838 template<typename _UIntType = unsigned int>
2839 class hypergeometric_distribution
2840 {
2841 static_assert(std::is_unsigned<_UIntType>::value, "template argument "
2842 "substituting _UIntType not an unsigned integral type");
2844 public:
2845 /** The type of the range of the distribution. */
2846 typedef _UIntType result_type;
2848 /** Parameter type. */
2849 struct param_type
2850 {
2851 typedef hypergeometric_distribution<_UIntType> distribution_type;
2852 friend class hypergeometric_distribution<_UIntType>;
2854 explicit
2855 param_type(result_type __N = 10, result_type __K = 5,
2856 result_type __n = 1)
2857 : _M_N{__N}, _M_K{__K}, _M_n{__n}
2858 {
2859 _GLIBCXX_DEBUG_ASSERT(_M_N >= _M_K);
2860 _GLIBCXX_DEBUG_ASSERT(_M_N >= _M_n);
2861 }
2863 result_type
2864 total_size() const
2865 { return _M_N; }
2867 result_type
2868 successful_size() const
2869 { return _M_K; }
2871 result_type
2872 unsuccessful_size() const
2873 { return _M_N - _M_K; }
2875 result_type
2876 total_draws() const
2877 { return _M_n; }
2879 friend bool
2880 operator==(const param_type& __p1, const param_type& __p2)
2881 { return (__p1._M_N == __p2._M_N)
2882 && (__p1._M_K == __p2._M_K)
2883 && (__p1._M_n == __p2._M_n); }
2885 private:
2887 result_type _M_N;
2888 result_type _M_K;
2889 result_type _M_n;
2890 };
2892 // constructors and member function
2893 explicit
2894 hypergeometric_distribution(result_type __N = 10, result_type __K = 5,
2895 result_type __n = 1)
2896 : _M_param{__N, __K, __n}
2897 { }
2899 explicit
2900 hypergeometric_distribution(const param_type& __p)
2901 : _M_param{__p}
2902 { }
2904 /**
2905 * @brief Resets the distribution state.
2906 */
2907 void
2908 reset()
2909 { }
2911 /**
2912 * @brief Returns the distribution parameter @p N,
2913 * the total number of items.
2914 */
2915 result_type
2916 total_size() const
2917 { return this->_M_param.total_size(); }
2919 /**
2920 * @brief Returns the distribution parameter @p K,
2921 * the total number of successful items.
2922 */
2923 result_type
2924 successful_size() const
2925 { return this->_M_param.successful_size(); }
2927 /**
2928 * @brief Returns the total number of unsuccessful items @f$ N - K @f$.
2929 */
2930 result_type
2931 unsuccessful_size() const
2932 { return this->_M_param.unsuccessful_size(); }
2934 /**
2935 * @brief Returns the distribution parameter @p n,
2936 * the total number of draws.
2937 */
2938 result_type
2939 total_draws() const
2940 { return this->_M_param.total_draws(); }
2942 /**
2943 * @brief Returns the parameter set of the distribution.
2944 */
2945 param_type
2946 param() const
2947 { return this->_M_param; }
2949 /**
2950 * @brief Sets the parameter set of the distribution.
2951 * @param __param The new parameter set of the distribution.
2952 */
2953 void
2954 param(const param_type& __param)
2955 { this->_M_param = __param; }
2957 /**
2958 * @brief Returns the greatest lower bound value of the distribution.
2959 */
2960 result_type
2961 min() const
2962 {
2963 using _IntType = typename std::make_signed<result_type>::type;
2964 return static_cast<result_type>(std::max(static_cast<_IntType>(0),
2965 static_cast<_IntType>(this->total_draws()
2966 - this->unsuccessful_size())));
2967 }
2969 /**
2970 * @brief Returns the least upper bound value of the distribution.
2971 */
2972 result_type
2973 max() const
2974 { return std::min(this->successful_size(), this->total_draws()); }
2976 /**
2977 * @brief Generating functions.
2978 */
2979 template<typename _UniformRandomNumberGenerator>
2980 result_type
2981 operator()(_UniformRandomNumberGenerator& __urng)
2982 { return this->operator()(__urng, this->_M_param); }
2984 template<typename _UniformRandomNumberGenerator>
2985 result_type
2986 operator()(_UniformRandomNumberGenerator& __urng,
2987 const param_type& __p);
2989 template<typename _ForwardIterator,
2990 typename _UniformRandomNumberGenerator>
2991 void
2992 __generate(_ForwardIterator __f, _ForwardIterator __t,
2993 _UniformRandomNumberGenerator& __urng)
2994 { this->__generate(__f, __t, __urng, this->_M_param); }
2996 template<typename _ForwardIterator,
2997 typename _UniformRandomNumberGenerator>
2998 void
2999 __generate(_ForwardIterator __f, _ForwardIterator __t,
3000 _UniformRandomNumberGenerator& __urng,
3001 const param_type& __p)
3002 { this->__generate_impl(__f, __t, __urng, __p); }
3004 template<typename _UniformRandomNumberGenerator>
3005 void
3006 __generate(result_type* __f, result_type* __t,
3007 _UniformRandomNumberGenerator& __urng,
3008 const param_type& __p)
3009 { this->__generate_impl(__f, __t, __urng, __p); }
3011 /**
3012 * @brief Return true if two hypergeometric distributions have the same
3013 * parameters and the sequences that would be generated
3014 * are equal.
3015 */
3016 friend bool
3017 operator==(const hypergeometric_distribution& __d1,
3018 const hypergeometric_distribution& __d2)
3019 { return __d1._M_param == __d2._M_param; }
3021 /**
3022 * @brief Inserts a %hypergeometric_distribution random number
3023 * distribution @p __x into the output stream @p __os.
3024 *
3025 * @param __os An output stream.
3026 * @param __x A %hypergeometric_distribution random number
3027 * distribution.
3028 *
3029 * @returns The output stream with the state of @p __x inserted or in
3030 * an error state.
3031 */
3032 template<typename _UIntType1, typename _CharT, typename _Traits>
3033 friend std::basic_ostream<_CharT, _Traits>&
3034 operator<<(std::basic_ostream<_CharT, _Traits>& __os,
3035 const __gnu_cxx::hypergeometric_distribution<_UIntType1>&
3036 __x);
3038 /**
3039 * @brief Extracts a %hypergeometric_distribution random number
3040 * distribution @p __x from the input stream @p __is.
3041 *
3042 * @param __is An input stream.
3043 * @param __x A %hypergeometric_distribution random number generator
3044 * distribution.
3045 *
3046 * @returns The input stream with @p __x extracted or in an error
3047 * state.
3048 */
3049 template<typename _UIntType1, typename _CharT, typename _Traits>
3050 friend std::basic_istream<_CharT, _Traits>&
3051 operator>>(std::basic_istream<_CharT, _Traits>& __is,
3052 __gnu_cxx::hypergeometric_distribution<_UIntType1>& __x);
3054 private:
3056 template<typename _ForwardIterator,
3057 typename _UniformRandomNumberGenerator>
3058 void
3059 __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
3060 _UniformRandomNumberGenerator& __urng,
3061 const param_type& __p);
3063 param_type _M_param;
3064 };
3066 /**
3067 * @brief Return true if two hypergeometric distributions are different.
3068 */
3069 template<typename _UIntType>
3070 inline bool
3071 operator!=(const __gnu_cxx::hypergeometric_distribution<_UIntType>& __d1,
3072 const __gnu_cxx::hypergeometric_distribution<_UIntType>& __d2)
3073 { return !(__d1 == __d2); }
3075 /**
3076 * @brief A logistic continuous distribution for random numbers.
3077 *
3078 * The formula for the logistic probability density function is
3079 * @f[
3080 * p(x|\a,\b) = \frac{e^{(x - a)/b}}{b[1 + e^{(x - a)/b}]^2}
3081 * @f]
3082 * where @f$b > 0@f$.
3083 *
3084 * The formula for the logistic probability function is
3085 * @f[
3086 * cdf(x|\a,\b) = \frac{e^{(x - a)/b}}{1 + e^{(x - a)/b}}
3087 * @f]
3088 * where @f$b > 0@f$.
3089 *
3090 * <table border=1 cellpadding=10 cellspacing=0>
3091 * <caption align=top>Distribution Statistics</caption>
3092 * <tr><td>Mean</td><td>@f$a@f$</td></tr>
3093 * <tr><td>Variance</td><td>@f$b^2\pi^2/3@f$</td></tr>
3094 * <tr><td>Range</td><td>@f$[0, \infty)@f$</td></tr>
3095 * </table>
3096 */
3097 template<typename _RealType = double>
3098 class
3099 logistic_distribution
3100 {
3101 static_assert(std::is_floating_point<_RealType>::value,
3102 "template argument not a floating point type");
3104 public:
3105 /** The type of the range of the distribution. */
3106 typedef _RealType result_type;
3107 /** Parameter type. */
3108 struct param_type
3109 {
3110 typedef logistic_distribution<result_type> distribution_type;
3112 param_type(result_type __a = result_type(0),
3113 result_type __b = result_type(1))
3114 : _M_a(__a), _M_b(__b)
3115 {
3116 _GLIBCXX_DEBUG_ASSERT(_M_b > result_type(0));
3117 }
3119 result_type
3120 a() const
3121 { return _M_a; }
3123 result_type
3124 b() const
3125 { return _M_b; }
3127 friend bool
3128 operator==(const param_type& __p1, const param_type& __p2)
3129 { return __p1._M_a == __p2._M_a
3130 && __p1._M_b == __p2._M_b; }
3132 private:
3133 void _M_initialize();
3135 result_type _M_a;
3136 result_type _M_b;
3137 };
3139 /**
3140 * @brief Constructors.
3141 */
3142 explicit
3143 logistic_distribution(result_type __a = result_type(0),
3144 result_type __b = result_type(1))
3145 : _M_param(__a, __b)
3146 { }
3148 explicit
3149 logistic_distribution(const param_type& __p)
3150 : _M_param(__p)
3151 { }
3153 /**
3154 * @brief Resets the distribution state.
3155 */
3156 void
3157 reset()
3158 { }
3160 /**
3161 * @brief Return the parameters of the distribution.
3162 */
3163 result_type
3164 a() const
3165 { return _M_param.a(); }
3167 result_type
3168 b() const
3169 { return _M_param.b(); }
3171 /**
3172 * @brief Returns the parameter set of the distribution.
3173 */
3174 param_type
3175 param() const
3176 { return _M_param; }
3178 /**
3179 * @brief Sets the parameter set of the distribution.
3180 * @param __param The new parameter set of the distribution.
3181 */
3182 void
3183 param(const param_type& __param)
3184 { _M_param = __param; }
3186 /**
3187 * @brief Returns the greatest lower bound value of the distribution.
3188 */
3189 result_type
3190 min() const
3191 { return -std::numeric_limits<result_type>::max(); }
3193 /**
3194 * @brief Returns the least upper bound value of the distribution.
3195 */
3196 result_type
3197 max() const
3198 { return std::numeric_limits<result_type>::max(); }
3200 /**
3201 * @brief Generating functions.
3202 */
3203 template<typename _UniformRandomNumberGenerator>
3204 result_type
3205 operator()(_UniformRandomNumberGenerator& __urng)
3206 { return this->operator()(__urng, this->_M_param); }
3208 template<typename _UniformRandomNumberGenerator>
3209 result_type
3210 operator()(_UniformRandomNumberGenerator&,
3211 const param_type&);
3213 template<typename _ForwardIterator,
3214 typename _UniformRandomNumberGenerator>
3215 void
3216 __generate(_ForwardIterator __f, _ForwardIterator __t,
3217 _UniformRandomNumberGenerator& __urng)
3218 { this->__generate(__f, __t, __urng, this->param()); }
3220 template<typename _ForwardIterator,
3221 typename _UniformRandomNumberGenerator>
3222 void
3223 __generate(_ForwardIterator __f, _ForwardIterator __t,
3224 _UniformRandomNumberGenerator& __urng,
3225 const param_type& __p)
3226 { this->__generate_impl(__f, __t, __urng, __p); }
3228 template<typename _UniformRandomNumberGenerator>
3229 void
3230 __generate(result_type* __f, result_type* __t,
3231 _UniformRandomNumberGenerator& __urng,
3232 const param_type& __p)
3233 { this->__generate_impl(__f, __t, __urng, __p); }
3235 /**
3236 * @brief Return true if two logistic distributions have
3237 * the same parameters and the sequences that would
3238 * be generated are equal.
3239 */
3240 template<typename _RealType1>
3241 friend bool
3242 operator==(const logistic_distribution<_RealType1>& __d1,
3243 const logistic_distribution<_RealType1>& __d2)
3244 { return __d1.param() == __d2.param(); }
3246 /**
3247 * @brief Inserts a %logistic_distribution random number distribution
3248 * @p __x into the output stream @p __os.
3249 *
3250 * @param __os An output stream.
3251 * @param __x A %logistic_distribution random number distribution.
3252 *
3253 * @returns The output stream with the state of @p __x inserted or in
3254 * an error state.
3255 */
3256 template<typename _RealType1, typename _CharT, typename _Traits>
3257 friend std::basic_ostream<_CharT, _Traits>&
3258 operator<<(std::basic_ostream<_CharT, _Traits>&,
3259 const logistic_distribution<_RealType1>&);
3261 /**
3262 * @brief Extracts a %logistic_distribution random number distribution
3263 * @p __x from the input stream @p __is.
3264 *
3265 * @param __is An input stream.
3266 * @param __x A %logistic_distribution random number
3267 * generator engine.
3268 *
3269 * @returns The input stream with @p __x extracted or in an error state.
3270 */
3271 template<typename _RealType1, typename _CharT, typename _Traits>
3272 friend std::basic_istream<_CharT, _Traits>&
3273 operator>>(std::basic_istream<_CharT, _Traits>&,
3274 logistic_distribution<_RealType1>&);
3276 private:
3277 template<typename _ForwardIterator,
3278 typename _UniformRandomNumberGenerator>
3279 void
3280 __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
3281 _UniformRandomNumberGenerator& __urng,
3282 const param_type& __p);
3284 param_type _M_param;
3285 };
3287 /**
3288 * @brief Return true if two logistic distributions are not equal.
3289 */
3290 template<typename _RealType1>
3291 inline bool
3292 operator!=(const logistic_distribution<_RealType1>& __d1,
3293 const logistic_distribution<_RealType1>& __d2)
3294 { return !(__d1 == __d2); }
3296 _GLIBCXX_END_NAMESPACE_VERSION
3297 } // namespace __gnu_cxx
3299 #include "ext/opt_random.h"
3300 #include "random.tcc"
3302 #endif // _GLIBCXX_USE_C99_STDINT_TR1
3304 #endif // C++11
3306 #endif // _EXT_RANDOM