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