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