xref: /netbsd-src/external/gpl3/gcc.old/dist/libstdc++-v3/include/ext/random.tcc (revision 82d56013d7b633d116a93943de88e08335357a7c)
1 // Random number extensions -*- C++ -*-
2 
3 // Copyright (C) 2012-2019 Free Software Foundation, Inc.
4 //
5 // This file is part of the GNU ISO C++ Library.  This library is free
6 // software; you can redistribute it and/or modify it under the
7 // terms of the GNU General Public License as published by the
8 // Free Software Foundation; either version 3, or (at your option)
9 // any later version.
10 
11 // This library is distributed in the hope that it will be useful,
12 // but WITHOUT ANY WARRANTY; without even the implied warranty of
13 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
14 // GNU General Public License for more details.
15 
16 // Under Section 7 of GPL version 3, you are granted additional
17 // permissions described in the GCC Runtime Library Exception, version
18 // 3.1, as published by the Free Software Foundation.
19 
20 // You should have received a copy of the GNU General Public License and
21 // a copy of the GCC Runtime Library Exception along with this program;
22 // see the files COPYING3 and COPYING.RUNTIME respectively.  If not, see
23 // <http://www.gnu.org/licenses/>.
24 
25 /** @file ext/random.tcc
26  *  This is an internal header file, included by other library headers.
27  *  Do not attempt to use it directly. @headername{ext/random}
28  */
29 
30 #ifndef _EXT_RANDOM_TCC
31 #define _EXT_RANDOM_TCC 1
32 
33 #pragma GCC system_header
34 
35 namespace __gnu_cxx _GLIBCXX_VISIBILITY(default)
36 {
37 _GLIBCXX_BEGIN_NAMESPACE_VERSION
38 
39 #if __BYTE_ORDER__ == __ORDER_LITTLE_ENDIAN__
40 
41   template<typename _UIntType, size_t __m,
42 	   size_t __pos1, size_t __sl1, size_t __sl2,
43 	   size_t __sr1, size_t __sr2,
44 	   uint32_t __msk1, uint32_t __msk2,
45 	   uint32_t __msk3, uint32_t __msk4,
46 	   uint32_t __parity1, uint32_t __parity2,
47 	   uint32_t __parity3, uint32_t __parity4>
48     void simd_fast_mersenne_twister_engine<_UIntType, __m,
49 					   __pos1, __sl1, __sl2, __sr1, __sr2,
50 					   __msk1, __msk2, __msk3, __msk4,
51 					   __parity1, __parity2, __parity3,
52 					   __parity4>::
53     seed(_UIntType __seed)
54     {
55       _M_state32[0] = static_cast<uint32_t>(__seed);
56       for (size_t __i = 1; __i < _M_nstate32; ++__i)
57 	_M_state32[__i] = (1812433253UL
58 			   * (_M_state32[__i - 1] ^ (_M_state32[__i - 1] >> 30))
59 			   + __i);
60       _M_pos = state_size;
61       _M_period_certification();
62     }
63 
64 
65   namespace {
66 
67     inline uint32_t _Func1(uint32_t __x)
68     {
69       return (__x ^ (__x >> 27)) * UINT32_C(1664525);
70     }
71 
72     inline uint32_t _Func2(uint32_t __x)
73     {
74       return (__x ^ (__x >> 27)) * UINT32_C(1566083941);
75     }
76 
77   }
78 
79 
80   template<typename _UIntType, size_t __m,
81 	   size_t __pos1, size_t __sl1, size_t __sl2,
82 	   size_t __sr1, size_t __sr2,
83 	   uint32_t __msk1, uint32_t __msk2,
84 	   uint32_t __msk3, uint32_t __msk4,
85 	   uint32_t __parity1, uint32_t __parity2,
86 	   uint32_t __parity3, uint32_t __parity4>
87     template<typename _Sseq>
88       auto
89       simd_fast_mersenne_twister_engine<_UIntType, __m,
90 					__pos1, __sl1, __sl2, __sr1, __sr2,
91 					__msk1, __msk2, __msk3, __msk4,
92 					__parity1, __parity2, __parity3,
93 					__parity4>::
94       seed(_Sseq& __q)
95       -> _If_seed_seq<_Sseq>
96       {
97 	size_t __lag;
98 
99 	if (_M_nstate32 >= 623)
100 	  __lag = 11;
101 	else if (_M_nstate32 >= 68)
102 	  __lag = 7;
103 	else if (_M_nstate32 >= 39)
104 	  __lag = 5;
105 	else
106 	  __lag = 3;
107 	const size_t __mid = (_M_nstate32 - __lag) / 2;
108 
109 	std::fill(_M_state32, _M_state32 + _M_nstate32, UINT32_C(0x8b8b8b8b));
110 	uint32_t __arr[_M_nstate32];
111 	__q.generate(__arr + 0, __arr + _M_nstate32);
112 
113 	uint32_t __r = _Func1(_M_state32[0] ^ _M_state32[__mid]
114 			      ^ _M_state32[_M_nstate32  - 1]);
115 	_M_state32[__mid] += __r;
116 	__r += _M_nstate32;
117 	_M_state32[__mid + __lag] += __r;
118 	_M_state32[0] = __r;
119 
120 	for (size_t __i = 1, __j = 0; __j < _M_nstate32; ++__j)
121 	  {
122 	    __r = _Func1(_M_state32[__i]
123 			 ^ _M_state32[(__i + __mid) % _M_nstate32]
124 			 ^ _M_state32[(__i + _M_nstate32 - 1) % _M_nstate32]);
125 	    _M_state32[(__i + __mid) % _M_nstate32] += __r;
126 	    __r += __arr[__j] + __i;
127 	    _M_state32[(__i + __mid + __lag) % _M_nstate32] += __r;
128 	    _M_state32[__i] = __r;
129 	    __i = (__i + 1) % _M_nstate32;
130 	  }
131 	for (size_t __j = 0; __j < _M_nstate32; ++__j)
132 	  {
133 	    const size_t __i = (__j + 1) % _M_nstate32;
134 	    __r = _Func2(_M_state32[__i]
135 			 + _M_state32[(__i + __mid) % _M_nstate32]
136 			 + _M_state32[(__i + _M_nstate32 - 1) % _M_nstate32]);
137 	    _M_state32[(__i + __mid) % _M_nstate32] ^= __r;
138 	    __r -= __i;
139 	    _M_state32[(__i + __mid + __lag) % _M_nstate32] ^= __r;
140 	    _M_state32[__i] = __r;
141 	  }
142 
143 	_M_pos = state_size;
144 	_M_period_certification();
145       }
146 
147 
148   template<typename _UIntType, size_t __m,
149 	   size_t __pos1, size_t __sl1, size_t __sl2,
150 	   size_t __sr1, size_t __sr2,
151 	   uint32_t __msk1, uint32_t __msk2,
152 	   uint32_t __msk3, uint32_t __msk4,
153 	   uint32_t __parity1, uint32_t __parity2,
154 	   uint32_t __parity3, uint32_t __parity4>
155     void simd_fast_mersenne_twister_engine<_UIntType, __m,
156 					   __pos1, __sl1, __sl2, __sr1, __sr2,
157 					   __msk1, __msk2, __msk3, __msk4,
158 					   __parity1, __parity2, __parity3,
159 					   __parity4>::
160     _M_period_certification(void)
161     {
162       static const uint32_t __parity[4] = { __parity1, __parity2,
163 					    __parity3, __parity4 };
164       uint32_t __inner = 0;
165       for (size_t __i = 0; __i < 4; ++__i)
166 	if (__parity[__i] != 0)
167 	  __inner ^= _M_state32[__i] & __parity[__i];
168 
169       if (__builtin_parity(__inner) & 1)
170 	return;
171       for (size_t __i = 0; __i < 4; ++__i)
172 	if (__parity[__i] != 0)
173 	  {
174 	    _M_state32[__i] ^= 1 << (__builtin_ffs(__parity[__i]) - 1);
175 	    return;
176 	  }
177       __builtin_unreachable();
178     }
179 
180 
181   template<typename _UIntType, size_t __m,
182 	   size_t __pos1, size_t __sl1, size_t __sl2,
183 	   size_t __sr1, size_t __sr2,
184 	   uint32_t __msk1, uint32_t __msk2,
185 	   uint32_t __msk3, uint32_t __msk4,
186 	   uint32_t __parity1, uint32_t __parity2,
187 	   uint32_t __parity3, uint32_t __parity4>
188     void simd_fast_mersenne_twister_engine<_UIntType, __m,
189 					   __pos1, __sl1, __sl2, __sr1, __sr2,
190 					   __msk1, __msk2, __msk3, __msk4,
191 					   __parity1, __parity2, __parity3,
192 					   __parity4>::
193     discard(unsigned long long __z)
194     {
195       while (__z > state_size - _M_pos)
196 	{
197 	  __z -= state_size - _M_pos;
198 
199 	  _M_gen_rand();
200 	}
201 
202       _M_pos += __z;
203     }
204 
205 
206 #ifndef  _GLIBCXX_OPT_HAVE_RANDOM_SFMT_GEN_READ
207 
208   namespace {
209 
210     template<size_t __shift>
211       inline void __rshift(uint32_t *__out, const uint32_t *__in)
212       {
213 	uint64_t __th = ((static_cast<uint64_t>(__in[3]) << 32)
214 			 | static_cast<uint64_t>(__in[2]));
215 	uint64_t __tl = ((static_cast<uint64_t>(__in[1]) << 32)
216 			 | static_cast<uint64_t>(__in[0]));
217 
218 	uint64_t __oh = __th >> (__shift * 8);
219 	uint64_t __ol = __tl >> (__shift * 8);
220 	__ol |= __th << (64 - __shift * 8);
221 	__out[1] = static_cast<uint32_t>(__ol >> 32);
222 	__out[0] = static_cast<uint32_t>(__ol);
223 	__out[3] = static_cast<uint32_t>(__oh >> 32);
224 	__out[2] = static_cast<uint32_t>(__oh);
225       }
226 
227 
228     template<size_t __shift>
229       inline void __lshift(uint32_t *__out, const uint32_t *__in)
230       {
231 	uint64_t __th = ((static_cast<uint64_t>(__in[3]) << 32)
232 			 | static_cast<uint64_t>(__in[2]));
233 	uint64_t __tl = ((static_cast<uint64_t>(__in[1]) << 32)
234 			 | static_cast<uint64_t>(__in[0]));
235 
236 	uint64_t __oh = __th << (__shift * 8);
237 	uint64_t __ol = __tl << (__shift * 8);
238 	__oh |= __tl >> (64 - __shift * 8);
239 	__out[1] = static_cast<uint32_t>(__ol >> 32);
240 	__out[0] = static_cast<uint32_t>(__ol);
241 	__out[3] = static_cast<uint32_t>(__oh >> 32);
242 	__out[2] = static_cast<uint32_t>(__oh);
243       }
244 
245 
246     template<size_t __sl1, size_t __sl2, size_t __sr1, size_t __sr2,
247 	     uint32_t __msk1, uint32_t __msk2, uint32_t __msk3, uint32_t __msk4>
248       inline void __recursion(uint32_t *__r,
249 			      const uint32_t *__a, const uint32_t *__b,
250 			      const uint32_t *__c, const uint32_t *__d)
251       {
252 	uint32_t __x[4];
253 	uint32_t __y[4];
254 
255 	__lshift<__sl2>(__x, __a);
256 	__rshift<__sr2>(__y, __c);
257 	__r[0] = (__a[0] ^ __x[0] ^ ((__b[0] >> __sr1) & __msk1)
258 		  ^ __y[0] ^ (__d[0] << __sl1));
259 	__r[1] = (__a[1] ^ __x[1] ^ ((__b[1] >> __sr1) & __msk2)
260 		  ^ __y[1] ^ (__d[1] << __sl1));
261 	__r[2] = (__a[2] ^ __x[2] ^ ((__b[2] >> __sr1) & __msk3)
262 		  ^ __y[2] ^ (__d[2] << __sl1));
263 	__r[3] = (__a[3] ^ __x[3] ^ ((__b[3] >> __sr1) & __msk4)
264 		  ^ __y[3] ^ (__d[3] << __sl1));
265       }
266 
267   }
268 
269 
270   template<typename _UIntType, size_t __m,
271 	   size_t __pos1, size_t __sl1, size_t __sl2,
272 	   size_t __sr1, size_t __sr2,
273 	   uint32_t __msk1, uint32_t __msk2,
274 	   uint32_t __msk3, uint32_t __msk4,
275 	   uint32_t __parity1, uint32_t __parity2,
276 	   uint32_t __parity3, uint32_t __parity4>
277     void simd_fast_mersenne_twister_engine<_UIntType, __m,
278 					   __pos1, __sl1, __sl2, __sr1, __sr2,
279 					   __msk1, __msk2, __msk3, __msk4,
280 					   __parity1, __parity2, __parity3,
281 					   __parity4>::
282     _M_gen_rand(void)
283     {
284       const uint32_t *__r1 = &_M_state32[_M_nstate32 - 8];
285       const uint32_t *__r2 = &_M_state32[_M_nstate32 - 4];
286       static constexpr size_t __pos1_32 = __pos1 * 4;
287 
288       size_t __i;
289       for (__i = 0; __i < _M_nstate32 - __pos1_32; __i += 4)
290 	{
291 	  __recursion<__sl1, __sl2, __sr1, __sr2,
292 		      __msk1, __msk2, __msk3, __msk4>
293 	    (&_M_state32[__i], &_M_state32[__i],
294 	     &_M_state32[__i + __pos1_32], __r1, __r2);
295 	  __r1 = __r2;
296 	  __r2 = &_M_state32[__i];
297 	}
298 
299       for (; __i < _M_nstate32; __i += 4)
300 	{
301 	  __recursion<__sl1, __sl2, __sr1, __sr2,
302 		      __msk1, __msk2, __msk3, __msk4>
303 	    (&_M_state32[__i], &_M_state32[__i],
304 	     &_M_state32[__i + __pos1_32 - _M_nstate32], __r1, __r2);
305 	  __r1 = __r2;
306 	  __r2 = &_M_state32[__i];
307 	}
308 
309       _M_pos = 0;
310     }
311 
312 #endif
313 
314 #ifndef _GLIBCXX_OPT_HAVE_RANDOM_SFMT_OPERATOREQUAL
315   template<typename _UIntType, size_t __m,
316 	   size_t __pos1, size_t __sl1, size_t __sl2,
317 	   size_t __sr1, size_t __sr2,
318 	   uint32_t __msk1, uint32_t __msk2,
319 	   uint32_t __msk3, uint32_t __msk4,
320 	   uint32_t __parity1, uint32_t __parity2,
321 	   uint32_t __parity3, uint32_t __parity4>
322     bool
323     operator==(const __gnu_cxx::simd_fast_mersenne_twister_engine<_UIntType,
324 	       __m, __pos1, __sl1, __sl2, __sr1, __sr2,
325 	       __msk1, __msk2, __msk3, __msk4,
326 	       __parity1, __parity2, __parity3, __parity4>& __lhs,
327 	       const __gnu_cxx::simd_fast_mersenne_twister_engine<_UIntType,
328 	       __m, __pos1, __sl1, __sl2, __sr1, __sr2,
329 	       __msk1, __msk2, __msk3, __msk4,
330 	       __parity1, __parity2, __parity3, __parity4>& __rhs)
331     {
332       typedef __gnu_cxx::simd_fast_mersenne_twister_engine<_UIntType,
333 	       __m, __pos1, __sl1, __sl2, __sr1, __sr2,
334 	       __msk1, __msk2, __msk3, __msk4,
335 	       __parity1, __parity2, __parity3, __parity4> __engine;
336       return (std::equal(__lhs._M_stateT,
337 			 __lhs._M_stateT + __engine::state_size,
338 			 __rhs._M_stateT)
339 	      && __lhs._M_pos == __rhs._M_pos);
340     }
341 #endif
342 
343   template<typename _UIntType, size_t __m,
344 	   size_t __pos1, size_t __sl1, size_t __sl2,
345 	   size_t __sr1, size_t __sr2,
346 	   uint32_t __msk1, uint32_t __msk2,
347 	   uint32_t __msk3, uint32_t __msk4,
348 	   uint32_t __parity1, uint32_t __parity2,
349 	   uint32_t __parity3, uint32_t __parity4,
350 	   typename _CharT, typename _Traits>
351     std::basic_ostream<_CharT, _Traits>&
352     operator<<(std::basic_ostream<_CharT, _Traits>& __os,
353 	       const __gnu_cxx::simd_fast_mersenne_twister_engine<_UIntType,
354 	       __m, __pos1, __sl1, __sl2, __sr1, __sr2,
355 	       __msk1, __msk2, __msk3, __msk4,
356 	       __parity1, __parity2, __parity3, __parity4>& __x)
357     {
358       typedef std::basic_ostream<_CharT, _Traits> __ostream_type;
359       typedef typename __ostream_type::ios_base __ios_base;
360 
361       const typename __ios_base::fmtflags __flags = __os.flags();
362       const _CharT __fill = __os.fill();
363       const _CharT __space = __os.widen(' ');
364       __os.flags(__ios_base::dec | __ios_base::fixed | __ios_base::left);
365       __os.fill(__space);
366 
367       for (size_t __i = 0; __i < __x._M_nstate32; ++__i)
368 	__os << __x._M_state32[__i] << __space;
369       __os << __x._M_pos;
370 
371       __os.flags(__flags);
372       __os.fill(__fill);
373       return __os;
374     }
375 
376 
377   template<typename _UIntType, size_t __m,
378 	   size_t __pos1, size_t __sl1, size_t __sl2,
379 	   size_t __sr1, size_t __sr2,
380 	   uint32_t __msk1, uint32_t __msk2,
381 	   uint32_t __msk3, uint32_t __msk4,
382 	   uint32_t __parity1, uint32_t __parity2,
383 	   uint32_t __parity3, uint32_t __parity4,
384 	   typename _CharT, typename _Traits>
385     std::basic_istream<_CharT, _Traits>&
386     operator>>(std::basic_istream<_CharT, _Traits>& __is,
387 	       __gnu_cxx::simd_fast_mersenne_twister_engine<_UIntType,
388 	       __m, __pos1, __sl1, __sl2, __sr1, __sr2,
389 	       __msk1, __msk2, __msk3, __msk4,
390 	       __parity1, __parity2, __parity3, __parity4>& __x)
391     {
392       typedef std::basic_istream<_CharT, _Traits> __istream_type;
393       typedef typename __istream_type::ios_base __ios_base;
394 
395       const typename __ios_base::fmtflags __flags = __is.flags();
396       __is.flags(__ios_base::dec | __ios_base::skipws);
397 
398       for (size_t __i = 0; __i < __x._M_nstate32; ++__i)
399 	__is >> __x._M_state32[__i];
400       __is >> __x._M_pos;
401 
402       __is.flags(__flags);
403       return __is;
404     }
405 
406 #endif // __BYTE_ORDER__ == __ORDER_LITTLE_ENDIAN__
407 
408   /**
409    * Iteration method due to M.D. J<o:>hnk.
410    *
411    * M.D. J<o:>hnk, Erzeugung von betaverteilten und gammaverteilten
412    * Zufallszahlen, Metrika, Volume 8, 1964
413    */
414   template<typename _RealType>
415     template<typename _UniformRandomNumberGenerator>
416       typename beta_distribution<_RealType>::result_type
417       beta_distribution<_RealType>::
418       operator()(_UniformRandomNumberGenerator& __urng,
419 		 const param_type& __param)
420       {
421 	std::__detail::_Adaptor<_UniformRandomNumberGenerator, result_type>
422 	  __aurng(__urng);
423 
424 	result_type __x, __y;
425 	do
426 	  {
427 	    __x = std::exp(std::log(__aurng()) / __param.alpha());
428 	    __y = std::exp(std::log(__aurng()) / __param.beta());
429 	  }
430 	while (__x + __y > result_type(1));
431 
432 	return __x / (__x + __y);
433       }
434 
435   template<typename _RealType>
436     template<typename _OutputIterator,
437 	     typename _UniformRandomNumberGenerator>
438       void
439       beta_distribution<_RealType>::
440       __generate_impl(_OutputIterator __f, _OutputIterator __t,
441 		      _UniformRandomNumberGenerator& __urng,
442 		      const param_type& __param)
443       {
444 	__glibcxx_function_requires(_OutputIteratorConcept<_OutputIterator,
445 	    result_type>)
446 
447 	std::__detail::_Adaptor<_UniformRandomNumberGenerator, result_type>
448 	  __aurng(__urng);
449 
450 	while (__f != __t)
451 	  {
452 	    result_type __x, __y;
453 	    do
454 	      {
455 		__x = std::exp(std::log(__aurng()) / __param.alpha());
456 		__y = std::exp(std::log(__aurng()) / __param.beta());
457 	      }
458 	    while (__x + __y > result_type(1));
459 
460 	    *__f++ = __x / (__x + __y);
461 	  }
462       }
463 
464   template<typename _RealType, typename _CharT, typename _Traits>
465     std::basic_ostream<_CharT, _Traits>&
466     operator<<(std::basic_ostream<_CharT, _Traits>& __os,
467 	       const __gnu_cxx::beta_distribution<_RealType>& __x)
468     {
469       typedef std::basic_ostream<_CharT, _Traits>  __ostream_type;
470       typedef typename __ostream_type::ios_base    __ios_base;
471 
472       const typename __ios_base::fmtflags __flags = __os.flags();
473       const _CharT __fill = __os.fill();
474       const std::streamsize __precision = __os.precision();
475       const _CharT __space = __os.widen(' ');
476       __os.flags(__ios_base::scientific | __ios_base::left);
477       __os.fill(__space);
478       __os.precision(std::numeric_limits<_RealType>::max_digits10);
479 
480       __os << __x.alpha() << __space << __x.beta();
481 
482       __os.flags(__flags);
483       __os.fill(__fill);
484       __os.precision(__precision);
485       return __os;
486     }
487 
488   template<typename _RealType, typename _CharT, typename _Traits>
489     std::basic_istream<_CharT, _Traits>&
490     operator>>(std::basic_istream<_CharT, _Traits>& __is,
491 	       __gnu_cxx::beta_distribution<_RealType>& __x)
492     {
493       typedef std::basic_istream<_CharT, _Traits>  __istream_type;
494       typedef typename __istream_type::ios_base    __ios_base;
495 
496       const typename __ios_base::fmtflags __flags = __is.flags();
497       __is.flags(__ios_base::dec | __ios_base::skipws);
498 
499       _RealType __alpha_val, __beta_val;
500       __is >> __alpha_val >> __beta_val;
501       __x.param(typename __gnu_cxx::beta_distribution<_RealType>::
502 		param_type(__alpha_val, __beta_val));
503 
504       __is.flags(__flags);
505       return __is;
506     }
507 
508 
509   template<std::size_t _Dimen, typename _RealType>
510     template<typename _InputIterator1, typename _InputIterator2>
511       void
512       normal_mv_distribution<_Dimen, _RealType>::param_type::
513       _M_init_full(_InputIterator1 __meanbegin, _InputIterator1 __meanend,
514 		   _InputIterator2 __varcovbegin, _InputIterator2 __varcovend)
515       {
516 	__glibcxx_function_requires(_InputIteratorConcept<_InputIterator1>)
517 	__glibcxx_function_requires(_InputIteratorConcept<_InputIterator2>)
518 	std::fill(std::copy(__meanbegin, __meanend, _M_mean.begin()),
519 		  _M_mean.end(), _RealType(0));
520 
521 	// Perform the Cholesky decomposition
522 	auto __w = _M_t.begin();
523 	for (size_t __j = 0; __j < _Dimen; ++__j)
524 	  {
525 	    _RealType __sum = _RealType(0);
526 
527 	    auto __slitbegin = __w;
528 	    auto __cit = _M_t.begin();
529 	    for (size_t __i = 0; __i < __j; ++__i)
530 	      {
531 		auto __slit = __slitbegin;
532 		_RealType __s = *__varcovbegin++;
533 		for (size_t __k = 0; __k < __i; ++__k)
534 		  __s -= *__slit++ * *__cit++;
535 
536 		*__w++ = __s /= *__cit++;
537 		__sum += __s * __s;
538 	      }
539 
540 	    __sum = *__varcovbegin - __sum;
541 	    if (__builtin_expect(__sum <= _RealType(0), 0))
542 	      std::__throw_runtime_error(__N("normal_mv_distribution::"
543 					     "param_type::_M_init_full"));
544 	    *__w++ = std::sqrt(__sum);
545 
546 	    std::advance(__varcovbegin, _Dimen - __j);
547 	  }
548       }
549 
550   template<std::size_t _Dimen, typename _RealType>
551     template<typename _InputIterator1, typename _InputIterator2>
552       void
553       normal_mv_distribution<_Dimen, _RealType>::param_type::
554       _M_init_lower(_InputIterator1 __meanbegin, _InputIterator1 __meanend,
555 		    _InputIterator2 __varcovbegin, _InputIterator2 __varcovend)
556       {
557 	__glibcxx_function_requires(_InputIteratorConcept<_InputIterator1>)
558 	__glibcxx_function_requires(_InputIteratorConcept<_InputIterator2>)
559 	std::fill(std::copy(__meanbegin, __meanend, _M_mean.begin()),
560 		  _M_mean.end(), _RealType(0));
561 
562 	// Perform the Cholesky decomposition
563 	auto __w = _M_t.begin();
564 	for (size_t __j = 0; __j < _Dimen; ++__j)
565 	  {
566 	    _RealType __sum = _RealType(0);
567 
568 	    auto __slitbegin = __w;
569 	    auto __cit = _M_t.begin();
570 	    for (size_t __i = 0; __i < __j; ++__i)
571 	      {
572 		auto __slit = __slitbegin;
573 		_RealType __s = *__varcovbegin++;
574 		for (size_t __k = 0; __k < __i; ++__k)
575 		  __s -= *__slit++ * *__cit++;
576 
577 		*__w++ = __s /= *__cit++;
578 		__sum += __s * __s;
579 	      }
580 
581 	    __sum = *__varcovbegin++ - __sum;
582 	    if (__builtin_expect(__sum <= _RealType(0), 0))
583 	      std::__throw_runtime_error(__N("normal_mv_distribution::"
584 					     "param_type::_M_init_full"));
585 	    *__w++ = std::sqrt(__sum);
586 	  }
587       }
588 
589   template<std::size_t _Dimen, typename _RealType>
590     template<typename _InputIterator1, typename _InputIterator2>
591       void
592       normal_mv_distribution<_Dimen, _RealType>::param_type::
593       _M_init_diagonal(_InputIterator1 __meanbegin, _InputIterator1 __meanend,
594 		       _InputIterator2 __varbegin, _InputIterator2 __varend)
595       {
596 	__glibcxx_function_requires(_InputIteratorConcept<_InputIterator1>)
597 	__glibcxx_function_requires(_InputIteratorConcept<_InputIterator2>)
598 	std::fill(std::copy(__meanbegin, __meanend, _M_mean.begin()),
599 		  _M_mean.end(), _RealType(0));
600 
601 	auto __w = _M_t.begin();
602 	size_t __step = 0;
603 	while (__varbegin != __varend)
604 	  {
605 	    std::fill_n(__w, __step, _RealType(0));
606 	    __w += __step++;
607 	    if (__builtin_expect(*__varbegin < _RealType(0), 0))
608 	      std::__throw_runtime_error(__N("normal_mv_distribution::"
609 					     "param_type::_M_init_diagonal"));
610 	    *__w++ = std::sqrt(*__varbegin++);
611 	  }
612       }
613 
614   template<std::size_t _Dimen, typename _RealType>
615     template<typename _UniformRandomNumberGenerator>
616       typename normal_mv_distribution<_Dimen, _RealType>::result_type
617       normal_mv_distribution<_Dimen, _RealType>::
618       operator()(_UniformRandomNumberGenerator& __urng,
619 		 const param_type& __param)
620       {
621 	result_type __ret;
622 
623 	_M_nd.__generate(__ret.begin(), __ret.end(), __urng);
624 
625 	auto __t_it = __param._M_t.crbegin();
626 	for (size_t __i = _Dimen; __i > 0; --__i)
627 	  {
628 	    _RealType __sum = _RealType(0);
629 	    for (size_t __j = __i; __j > 0; --__j)
630 	      __sum += __ret[__j - 1] * *__t_it++;
631 	    __ret[__i - 1] = __sum;
632 	  }
633 
634 	return __ret;
635       }
636 
637   template<std::size_t _Dimen, typename _RealType>
638     template<typename _ForwardIterator, typename _UniformRandomNumberGenerator>
639       void
640       normal_mv_distribution<_Dimen, _RealType>::
641       __generate_impl(_ForwardIterator __f, _ForwardIterator __t,
642 		      _UniformRandomNumberGenerator& __urng,
643 		      const param_type& __param)
644       {
645 	__glibcxx_function_requires(_Mutable_ForwardIteratorConcept<
646 				    _ForwardIterator>)
647 	while (__f != __t)
648 	  *__f++ = this->operator()(__urng, __param);
649       }
650 
651   template<size_t _Dimen, typename _RealType>
652     bool
653     operator==(const __gnu_cxx::normal_mv_distribution<_Dimen, _RealType>&
654 	       __d1,
655 	       const __gnu_cxx::normal_mv_distribution<_Dimen, _RealType>&
656 	       __d2)
657     {
658       return __d1._M_param == __d2._M_param && __d1._M_nd == __d2._M_nd;
659     }
660 
661   template<size_t _Dimen, typename _RealType, typename _CharT, typename _Traits>
662     std::basic_ostream<_CharT, _Traits>&
663     operator<<(std::basic_ostream<_CharT, _Traits>& __os,
664 	       const __gnu_cxx::normal_mv_distribution<_Dimen, _RealType>& __x)
665     {
666       typedef std::basic_ostream<_CharT, _Traits>  __ostream_type;
667       typedef typename __ostream_type::ios_base    __ios_base;
668 
669       const typename __ios_base::fmtflags __flags = __os.flags();
670       const _CharT __fill = __os.fill();
671       const std::streamsize __precision = __os.precision();
672       const _CharT __space = __os.widen(' ');
673       __os.flags(__ios_base::scientific | __ios_base::left);
674       __os.fill(__space);
675       __os.precision(std::numeric_limits<_RealType>::max_digits10);
676 
677       auto __mean = __x._M_param.mean();
678       for (auto __it : __mean)
679 	__os << __it << __space;
680       auto __t = __x._M_param.varcov();
681       for (auto __it : __t)
682 	__os << __it << __space;
683 
684       __os << __x._M_nd;
685 
686       __os.flags(__flags);
687       __os.fill(__fill);
688       __os.precision(__precision);
689       return __os;
690     }
691 
692   template<size_t _Dimen, typename _RealType, typename _CharT, typename _Traits>
693     std::basic_istream<_CharT, _Traits>&
694     operator>>(std::basic_istream<_CharT, _Traits>& __is,
695 	       __gnu_cxx::normal_mv_distribution<_Dimen, _RealType>& __x)
696     {
697       typedef std::basic_istream<_CharT, _Traits>  __istream_type;
698       typedef typename __istream_type::ios_base    __ios_base;
699 
700       const typename __ios_base::fmtflags __flags = __is.flags();
701       __is.flags(__ios_base::dec | __ios_base::skipws);
702 
703       std::array<_RealType, _Dimen> __mean;
704       for (auto& __it : __mean)
705 	__is >> __it;
706       std::array<_RealType, _Dimen * (_Dimen + 1) / 2> __varcov;
707       for (auto& __it : __varcov)
708 	__is >> __it;
709 
710       __is >> __x._M_nd;
711 
712       __x.param(typename normal_mv_distribution<_Dimen, _RealType>::
713 		param_type(__mean.begin(), __mean.end(),
714 			   __varcov.begin(), __varcov.end()));
715 
716       __is.flags(__flags);
717       return __is;
718     }
719 
720 
721   template<typename _RealType>
722     template<typename _OutputIterator,
723 	     typename _UniformRandomNumberGenerator>
724       void
725       rice_distribution<_RealType>::
726       __generate_impl(_OutputIterator __f, _OutputIterator __t,
727 		      _UniformRandomNumberGenerator& __urng,
728 		      const param_type& __p)
729       {
730 	__glibcxx_function_requires(_OutputIteratorConcept<_OutputIterator,
731 	    result_type>)
732 
733 	while (__f != __t)
734 	  {
735 	    typename std::normal_distribution<result_type>::param_type
736 	      __px(__p.nu(), __p.sigma()), __py(result_type(0), __p.sigma());
737 	    result_type __x = this->_M_ndx(__px, __urng);
738 	    result_type __y = this->_M_ndy(__py, __urng);
739 #if _GLIBCXX_USE_C99_MATH_TR1
740 	    *__f++ = std::hypot(__x, __y);
741 #else
742 	    *__f++ = std::sqrt(__x * __x + __y * __y);
743 #endif
744 	  }
745       }
746 
747   template<typename _RealType, typename _CharT, typename _Traits>
748     std::basic_ostream<_CharT, _Traits>&
749     operator<<(std::basic_ostream<_CharT, _Traits>& __os,
750 	       const rice_distribution<_RealType>& __x)
751     {
752       typedef std::basic_ostream<_CharT, _Traits>  __ostream_type;
753       typedef typename __ostream_type::ios_base    __ios_base;
754 
755       const typename __ios_base::fmtflags __flags = __os.flags();
756       const _CharT __fill = __os.fill();
757       const std::streamsize __precision = __os.precision();
758       const _CharT __space = __os.widen(' ');
759       __os.flags(__ios_base::scientific | __ios_base::left);
760       __os.fill(__space);
761       __os.precision(std::numeric_limits<_RealType>::max_digits10);
762 
763       __os << __x.nu() << __space << __x.sigma();
764       __os << __space << __x._M_ndx;
765       __os << __space << __x._M_ndy;
766 
767       __os.flags(__flags);
768       __os.fill(__fill);
769       __os.precision(__precision);
770       return __os;
771     }
772 
773   template<typename _RealType, typename _CharT, typename _Traits>
774     std::basic_istream<_CharT, _Traits>&
775     operator>>(std::basic_istream<_CharT, _Traits>& __is,
776 	       rice_distribution<_RealType>& __x)
777     {
778       typedef std::basic_istream<_CharT, _Traits>  __istream_type;
779       typedef typename __istream_type::ios_base    __ios_base;
780 
781       const typename __ios_base::fmtflags __flags = __is.flags();
782       __is.flags(__ios_base::dec | __ios_base::skipws);
783 
784       _RealType __nu_val, __sigma_val;
785       __is >> __nu_val >> __sigma_val;
786       __is >> __x._M_ndx;
787       __is >> __x._M_ndy;
788       __x.param(typename rice_distribution<_RealType>::
789 		param_type(__nu_val, __sigma_val));
790 
791       __is.flags(__flags);
792       return __is;
793     }
794 
795 
796   template<typename _RealType>
797     template<typename _OutputIterator,
798 	     typename _UniformRandomNumberGenerator>
799       void
800       nakagami_distribution<_RealType>::
801       __generate_impl(_OutputIterator __f, _OutputIterator __t,
802 		      _UniformRandomNumberGenerator& __urng,
803 		      const param_type& __p)
804       {
805 	__glibcxx_function_requires(_OutputIteratorConcept<_OutputIterator,
806 	    result_type>)
807 
808 	typename std::gamma_distribution<result_type>::param_type
809 	  __pg(__p.mu(), __p.omega() / __p.mu());
810 	while (__f != __t)
811 	  *__f++ = std::sqrt(this->_M_gd(__pg, __urng));
812       }
813 
814   template<typename _RealType, typename _CharT, typename _Traits>
815     std::basic_ostream<_CharT, _Traits>&
816     operator<<(std::basic_ostream<_CharT, _Traits>& __os,
817 	       const nakagami_distribution<_RealType>& __x)
818     {
819       typedef std::basic_ostream<_CharT, _Traits>  __ostream_type;
820       typedef typename __ostream_type::ios_base    __ios_base;
821 
822       const typename __ios_base::fmtflags __flags = __os.flags();
823       const _CharT __fill = __os.fill();
824       const std::streamsize __precision = __os.precision();
825       const _CharT __space = __os.widen(' ');
826       __os.flags(__ios_base::scientific | __ios_base::left);
827       __os.fill(__space);
828       __os.precision(std::numeric_limits<_RealType>::max_digits10);
829 
830       __os << __x.mu() << __space << __x.omega();
831       __os << __space << __x._M_gd;
832 
833       __os.flags(__flags);
834       __os.fill(__fill);
835       __os.precision(__precision);
836       return __os;
837     }
838 
839   template<typename _RealType, typename _CharT, typename _Traits>
840     std::basic_istream<_CharT, _Traits>&
841     operator>>(std::basic_istream<_CharT, _Traits>& __is,
842 	       nakagami_distribution<_RealType>& __x)
843     {
844       typedef std::basic_istream<_CharT, _Traits>  __istream_type;
845       typedef typename __istream_type::ios_base    __ios_base;
846 
847       const typename __ios_base::fmtflags __flags = __is.flags();
848       __is.flags(__ios_base::dec | __ios_base::skipws);
849 
850       _RealType __mu_val, __omega_val;
851       __is >> __mu_val >> __omega_val;
852       __is >> __x._M_gd;
853       __x.param(typename nakagami_distribution<_RealType>::
854 		param_type(__mu_val, __omega_val));
855 
856       __is.flags(__flags);
857       return __is;
858     }
859 
860 
861   template<typename _RealType>
862     template<typename _OutputIterator,
863 	     typename _UniformRandomNumberGenerator>
864       void
865       pareto_distribution<_RealType>::
866       __generate_impl(_OutputIterator __f, _OutputIterator __t,
867 		      _UniformRandomNumberGenerator& __urng,
868 		      const param_type& __p)
869       {
870 	__glibcxx_function_requires(_OutputIteratorConcept<_OutputIterator,
871 	    result_type>)
872 
873 	result_type __mu_val = __p.mu();
874 	result_type __malphinv = -result_type(1) / __p.alpha();
875 	while (__f != __t)
876 	  *__f++ = __mu_val * std::pow(this->_M_ud(__urng), __malphinv);
877       }
878 
879   template<typename _RealType, typename _CharT, typename _Traits>
880     std::basic_ostream<_CharT, _Traits>&
881     operator<<(std::basic_ostream<_CharT, _Traits>& __os,
882 	       const pareto_distribution<_RealType>& __x)
883     {
884       typedef std::basic_ostream<_CharT, _Traits>  __ostream_type;
885       typedef typename __ostream_type::ios_base    __ios_base;
886 
887       const typename __ios_base::fmtflags __flags = __os.flags();
888       const _CharT __fill = __os.fill();
889       const std::streamsize __precision = __os.precision();
890       const _CharT __space = __os.widen(' ');
891       __os.flags(__ios_base::scientific | __ios_base::left);
892       __os.fill(__space);
893       __os.precision(std::numeric_limits<_RealType>::max_digits10);
894 
895       __os << __x.alpha() << __space << __x.mu();
896       __os << __space << __x._M_ud;
897 
898       __os.flags(__flags);
899       __os.fill(__fill);
900       __os.precision(__precision);
901       return __os;
902     }
903 
904   template<typename _RealType, typename _CharT, typename _Traits>
905     std::basic_istream<_CharT, _Traits>&
906     operator>>(std::basic_istream<_CharT, _Traits>& __is,
907 	       pareto_distribution<_RealType>& __x)
908     {
909       typedef std::basic_istream<_CharT, _Traits>  __istream_type;
910       typedef typename __istream_type::ios_base    __ios_base;
911 
912       const typename __ios_base::fmtflags __flags = __is.flags();
913       __is.flags(__ios_base::dec | __ios_base::skipws);
914 
915       _RealType __alpha_val, __mu_val;
916       __is >> __alpha_val >> __mu_val;
917       __is >> __x._M_ud;
918       __x.param(typename pareto_distribution<_RealType>::
919 		param_type(__alpha_val, __mu_val));
920 
921       __is.flags(__flags);
922       return __is;
923     }
924 
925 
926   template<typename _RealType>
927     template<typename _UniformRandomNumberGenerator>
928       typename k_distribution<_RealType>::result_type
929       k_distribution<_RealType>::
930       operator()(_UniformRandomNumberGenerator& __urng)
931       {
932 	result_type __x = this->_M_gd1(__urng);
933 	result_type __y = this->_M_gd2(__urng);
934 	return std::sqrt(__x * __y);
935       }
936 
937   template<typename _RealType>
938     template<typename _UniformRandomNumberGenerator>
939       typename k_distribution<_RealType>::result_type
940       k_distribution<_RealType>::
941       operator()(_UniformRandomNumberGenerator& __urng,
942 		 const param_type& __p)
943       {
944 	typename std::gamma_distribution<result_type>::param_type
945 	  __p1(__p.lambda(), result_type(1) / __p.lambda()),
946 	  __p2(__p.nu(), __p.mu() / __p.nu());
947 	result_type __x = this->_M_gd1(__p1, __urng);
948 	result_type __y = this->_M_gd2(__p2, __urng);
949 	return std::sqrt(__x * __y);
950       }
951 
952   template<typename _RealType>
953     template<typename _OutputIterator,
954 	     typename _UniformRandomNumberGenerator>
955       void
956       k_distribution<_RealType>::
957       __generate_impl(_OutputIterator __f, _OutputIterator __t,
958 		      _UniformRandomNumberGenerator& __urng,
959 		      const param_type& __p)
960       {
961 	__glibcxx_function_requires(_OutputIteratorConcept<_OutputIterator,
962 	    result_type>)
963 
964 	typename std::gamma_distribution<result_type>::param_type
965 	  __p1(__p.lambda(), result_type(1) / __p.lambda()),
966 	  __p2(__p.nu(), __p.mu() / __p.nu());
967 	while (__f != __t)
968 	  {
969 	    result_type __x = this->_M_gd1(__p1, __urng);
970 	    result_type __y = this->_M_gd2(__p2, __urng);
971 	    *__f++ = std::sqrt(__x * __y);
972 	  }
973       }
974 
975   template<typename _RealType, typename _CharT, typename _Traits>
976     std::basic_ostream<_CharT, _Traits>&
977     operator<<(std::basic_ostream<_CharT, _Traits>& __os,
978 	       const k_distribution<_RealType>& __x)
979     {
980       typedef std::basic_ostream<_CharT, _Traits>  __ostream_type;
981       typedef typename __ostream_type::ios_base    __ios_base;
982 
983       const typename __ios_base::fmtflags __flags = __os.flags();
984       const _CharT __fill = __os.fill();
985       const std::streamsize __precision = __os.precision();
986       const _CharT __space = __os.widen(' ');
987       __os.flags(__ios_base::scientific | __ios_base::left);
988       __os.fill(__space);
989       __os.precision(std::numeric_limits<_RealType>::max_digits10);
990 
991       __os << __x.lambda() << __space << __x.mu() << __space << __x.nu();
992       __os << __space << __x._M_gd1;
993       __os << __space << __x._M_gd2;
994 
995       __os.flags(__flags);
996       __os.fill(__fill);
997       __os.precision(__precision);
998       return __os;
999     }
1000 
1001   template<typename _RealType, typename _CharT, typename _Traits>
1002     std::basic_istream<_CharT, _Traits>&
1003     operator>>(std::basic_istream<_CharT, _Traits>& __is,
1004 	       k_distribution<_RealType>& __x)
1005     {
1006       typedef std::basic_istream<_CharT, _Traits>  __istream_type;
1007       typedef typename __istream_type::ios_base    __ios_base;
1008 
1009       const typename __ios_base::fmtflags __flags = __is.flags();
1010       __is.flags(__ios_base::dec | __ios_base::skipws);
1011 
1012       _RealType __lambda_val, __mu_val, __nu_val;
1013       __is >> __lambda_val >> __mu_val >> __nu_val;
1014       __is >> __x._M_gd1;
1015       __is >> __x._M_gd2;
1016       __x.param(typename k_distribution<_RealType>::
1017 		param_type(__lambda_val, __mu_val, __nu_val));
1018 
1019       __is.flags(__flags);
1020       return __is;
1021     }
1022 
1023 
1024   template<typename _RealType>
1025     template<typename _OutputIterator,
1026 	     typename _UniformRandomNumberGenerator>
1027       void
1028       arcsine_distribution<_RealType>::
1029       __generate_impl(_OutputIterator __f, _OutputIterator __t,
1030 		      _UniformRandomNumberGenerator& __urng,
1031 		      const param_type& __p)
1032       {
1033 	__glibcxx_function_requires(_OutputIteratorConcept<_OutputIterator,
1034 	    result_type>)
1035 
1036 	result_type __dif = __p.b() - __p.a();
1037 	result_type __sum = __p.a() + __p.b();
1038 	while (__f != __t)
1039 	  {
1040 	    result_type __x = std::sin(this->_M_ud(__urng));
1041 	    *__f++ = (__x * __dif + __sum) / result_type(2);
1042 	  }
1043       }
1044 
1045   template<typename _RealType, typename _CharT, typename _Traits>
1046     std::basic_ostream<_CharT, _Traits>&
1047     operator<<(std::basic_ostream<_CharT, _Traits>& __os,
1048 	       const arcsine_distribution<_RealType>& __x)
1049     {
1050       typedef std::basic_ostream<_CharT, _Traits>  __ostream_type;
1051       typedef typename __ostream_type::ios_base    __ios_base;
1052 
1053       const typename __ios_base::fmtflags __flags = __os.flags();
1054       const _CharT __fill = __os.fill();
1055       const std::streamsize __precision = __os.precision();
1056       const _CharT __space = __os.widen(' ');
1057       __os.flags(__ios_base::scientific | __ios_base::left);
1058       __os.fill(__space);
1059       __os.precision(std::numeric_limits<_RealType>::max_digits10);
1060 
1061       __os << __x.a() << __space << __x.b();
1062       __os << __space << __x._M_ud;
1063 
1064       __os.flags(__flags);
1065       __os.fill(__fill);
1066       __os.precision(__precision);
1067       return __os;
1068     }
1069 
1070   template<typename _RealType, typename _CharT, typename _Traits>
1071     std::basic_istream<_CharT, _Traits>&
1072     operator>>(std::basic_istream<_CharT, _Traits>& __is,
1073 	       arcsine_distribution<_RealType>& __x)
1074     {
1075       typedef std::basic_istream<_CharT, _Traits>  __istream_type;
1076       typedef typename __istream_type::ios_base    __ios_base;
1077 
1078       const typename __ios_base::fmtflags __flags = __is.flags();
1079       __is.flags(__ios_base::dec | __ios_base::skipws);
1080 
1081       _RealType __a, __b;
1082       __is >> __a >> __b;
1083       __is >> __x._M_ud;
1084       __x.param(typename arcsine_distribution<_RealType>::
1085 		param_type(__a, __b));
1086 
1087       __is.flags(__flags);
1088       return __is;
1089     }
1090 
1091 
1092   template<typename _RealType>
1093     template<typename _UniformRandomNumberGenerator>
1094       typename hoyt_distribution<_RealType>::result_type
1095       hoyt_distribution<_RealType>::
1096       operator()(_UniformRandomNumberGenerator& __urng)
1097       {
1098 	result_type __x = this->_M_ad(__urng);
1099 	result_type __y = this->_M_ed(__urng);
1100 	return (result_type(2) * this->q()
1101 		  / (result_type(1) + this->q() * this->q()))
1102 	       * std::sqrt(this->omega() * __x * __y);
1103       }
1104 
1105   template<typename _RealType>
1106     template<typename _UniformRandomNumberGenerator>
1107       typename hoyt_distribution<_RealType>::result_type
1108       hoyt_distribution<_RealType>::
1109       operator()(_UniformRandomNumberGenerator& __urng,
1110 		 const param_type& __p)
1111       {
1112 	result_type __q2 = __p.q() * __p.q();
1113 	result_type __num = result_type(0.5L) * (result_type(1) + __q2);
1114 	typename __gnu_cxx::arcsine_distribution<result_type>::param_type
1115 	  __pa(__num, __num / __q2);
1116 	result_type __x = this->_M_ad(__pa, __urng);
1117 	result_type __y = this->_M_ed(__urng);
1118 	return (result_type(2) * __p.q() / (result_type(1) + __q2))
1119 	       * std::sqrt(__p.omega() * __x * __y);
1120       }
1121 
1122   template<typename _RealType>
1123     template<typename _OutputIterator,
1124 	     typename _UniformRandomNumberGenerator>
1125       void
1126       hoyt_distribution<_RealType>::
1127       __generate_impl(_OutputIterator __f, _OutputIterator __t,
1128 		      _UniformRandomNumberGenerator& __urng,
1129 		      const param_type& __p)
1130       {
1131 	__glibcxx_function_requires(_OutputIteratorConcept<_OutputIterator,
1132 	    result_type>)
1133 
1134 	result_type __2q = result_type(2) * __p.q();
1135 	result_type __q2 = __p.q() * __p.q();
1136 	result_type __q2p1 = result_type(1) + __q2;
1137 	result_type __num = result_type(0.5L) * __q2p1;
1138 	result_type __omega = __p.omega();
1139 	typename __gnu_cxx::arcsine_distribution<result_type>::param_type
1140 	  __pa(__num, __num / __q2);
1141 	while (__f != __t)
1142 	  {
1143 	    result_type __x = this->_M_ad(__pa, __urng);
1144 	    result_type __y = this->_M_ed(__urng);
1145 	    *__f++ = (__2q / __q2p1) * std::sqrt(__omega * __x * __y);
1146 	  }
1147       }
1148 
1149   template<typename _RealType, typename _CharT, typename _Traits>
1150     std::basic_ostream<_CharT, _Traits>&
1151     operator<<(std::basic_ostream<_CharT, _Traits>& __os,
1152 	       const hoyt_distribution<_RealType>& __x)
1153     {
1154       typedef std::basic_ostream<_CharT, _Traits>  __ostream_type;
1155       typedef typename __ostream_type::ios_base    __ios_base;
1156 
1157       const typename __ios_base::fmtflags __flags = __os.flags();
1158       const _CharT __fill = __os.fill();
1159       const std::streamsize __precision = __os.precision();
1160       const _CharT __space = __os.widen(' ');
1161       __os.flags(__ios_base::scientific | __ios_base::left);
1162       __os.fill(__space);
1163       __os.precision(std::numeric_limits<_RealType>::max_digits10);
1164 
1165       __os << __x.q() << __space << __x.omega();
1166       __os << __space << __x._M_ad;
1167       __os << __space << __x._M_ed;
1168 
1169       __os.flags(__flags);
1170       __os.fill(__fill);
1171       __os.precision(__precision);
1172       return __os;
1173     }
1174 
1175   template<typename _RealType, typename _CharT, typename _Traits>
1176     std::basic_istream<_CharT, _Traits>&
1177     operator>>(std::basic_istream<_CharT, _Traits>& __is,
1178 	       hoyt_distribution<_RealType>& __x)
1179     {
1180       typedef std::basic_istream<_CharT, _Traits>  __istream_type;
1181       typedef typename __istream_type::ios_base    __ios_base;
1182 
1183       const typename __ios_base::fmtflags __flags = __is.flags();
1184       __is.flags(__ios_base::dec | __ios_base::skipws);
1185 
1186       _RealType __q, __omega;
1187       __is >> __q >> __omega;
1188       __is >> __x._M_ad;
1189       __is >> __x._M_ed;
1190       __x.param(typename hoyt_distribution<_RealType>::
1191 		param_type(__q, __omega));
1192 
1193       __is.flags(__flags);
1194       return __is;
1195     }
1196 
1197 
1198   template<typename _RealType>
1199     template<typename _OutputIterator,
1200 	     typename _UniformRandomNumberGenerator>
1201       void
1202       triangular_distribution<_RealType>::
1203       __generate_impl(_OutputIterator __f, _OutputIterator __t,
1204 		      _UniformRandomNumberGenerator& __urng,
1205 		      const param_type& __param)
1206       {
1207 	__glibcxx_function_requires(_OutputIteratorConcept<_OutputIterator,
1208 	    result_type>)
1209 
1210 	while (__f != __t)
1211 	  *__f++ = this->operator()(__urng, __param);
1212       }
1213 
1214   template<typename _RealType, typename _CharT, typename _Traits>
1215     std::basic_ostream<_CharT, _Traits>&
1216     operator<<(std::basic_ostream<_CharT, _Traits>& __os,
1217 	       const __gnu_cxx::triangular_distribution<_RealType>& __x)
1218     {
1219       typedef std::basic_ostream<_CharT, _Traits>  __ostream_type;
1220       typedef typename __ostream_type::ios_base    __ios_base;
1221 
1222       const typename __ios_base::fmtflags __flags = __os.flags();
1223       const _CharT __fill = __os.fill();
1224       const std::streamsize __precision = __os.precision();
1225       const _CharT __space = __os.widen(' ');
1226       __os.flags(__ios_base::scientific | __ios_base::left);
1227       __os.fill(__space);
1228       __os.precision(std::numeric_limits<_RealType>::max_digits10);
1229 
1230       __os << __x.a() << __space << __x.b() << __space << __x.c();
1231 
1232       __os.flags(__flags);
1233       __os.fill(__fill);
1234       __os.precision(__precision);
1235       return __os;
1236     }
1237 
1238   template<typename _RealType, typename _CharT, typename _Traits>
1239     std::basic_istream<_CharT, _Traits>&
1240     operator>>(std::basic_istream<_CharT, _Traits>& __is,
1241 	       __gnu_cxx::triangular_distribution<_RealType>& __x)
1242     {
1243       typedef std::basic_istream<_CharT, _Traits>  __istream_type;
1244       typedef typename __istream_type::ios_base    __ios_base;
1245 
1246       const typename __ios_base::fmtflags __flags = __is.flags();
1247       __is.flags(__ios_base::dec | __ios_base::skipws);
1248 
1249       _RealType __a, __b, __c;
1250       __is >> __a >> __b >> __c;
1251       __x.param(typename __gnu_cxx::triangular_distribution<_RealType>::
1252 		param_type(__a, __b, __c));
1253 
1254       __is.flags(__flags);
1255       return __is;
1256     }
1257 
1258 
1259   template<typename _RealType>
1260     template<typename _UniformRandomNumberGenerator>
1261       typename von_mises_distribution<_RealType>::result_type
1262       von_mises_distribution<_RealType>::
1263       operator()(_UniformRandomNumberGenerator& __urng,
1264 		 const param_type& __p)
1265       {
1266 	const result_type __pi
1267 	  = __gnu_cxx::__math_constants<result_type>::__pi;
1268 	std::__detail::_Adaptor<_UniformRandomNumberGenerator, result_type>
1269 	  __aurng(__urng);
1270 
1271 	result_type __f;
1272 	while (1)
1273 	  {
1274 	    result_type __rnd = std::cos(__pi * __aurng());
1275 	    __f = (result_type(1) + __p._M_r * __rnd) / (__p._M_r + __rnd);
1276 	    result_type __c = __p._M_kappa * (__p._M_r - __f);
1277 
1278 	    result_type __rnd2 = __aurng();
1279 	    if (__c * (result_type(2) - __c) > __rnd2)
1280 	      break;
1281 	    if (std::log(__c / __rnd2) >= __c - result_type(1))
1282 	      break;
1283 	  }
1284 
1285 	result_type __res = std::acos(__f);
1286 #if _GLIBCXX_USE_C99_MATH_TR1
1287 	__res = std::copysign(__res, __aurng() - result_type(0.5));
1288 #else
1289 	if (__aurng() < result_type(0.5))
1290 	  __res = -__res;
1291 #endif
1292 	__res += __p._M_mu;
1293 	if (__res > __pi)
1294 	  __res -= result_type(2) * __pi;
1295 	else if (__res < -__pi)
1296 	  __res += result_type(2) * __pi;
1297 	return __res;
1298       }
1299 
1300   template<typename _RealType>
1301     template<typename _OutputIterator,
1302 	     typename _UniformRandomNumberGenerator>
1303       void
1304       von_mises_distribution<_RealType>::
1305       __generate_impl(_OutputIterator __f, _OutputIterator __t,
1306 		      _UniformRandomNumberGenerator& __urng,
1307 		      const param_type& __param)
1308       {
1309 	__glibcxx_function_requires(_OutputIteratorConcept<_OutputIterator,
1310 	    result_type>)
1311 
1312 	while (__f != __t)
1313 	  *__f++ = this->operator()(__urng, __param);
1314       }
1315 
1316   template<typename _RealType, typename _CharT, typename _Traits>
1317     std::basic_ostream<_CharT, _Traits>&
1318     operator<<(std::basic_ostream<_CharT, _Traits>& __os,
1319 	       const __gnu_cxx::von_mises_distribution<_RealType>& __x)
1320     {
1321       typedef std::basic_ostream<_CharT, _Traits>  __ostream_type;
1322       typedef typename __ostream_type::ios_base    __ios_base;
1323 
1324       const typename __ios_base::fmtflags __flags = __os.flags();
1325       const _CharT __fill = __os.fill();
1326       const std::streamsize __precision = __os.precision();
1327       const _CharT __space = __os.widen(' ');
1328       __os.flags(__ios_base::scientific | __ios_base::left);
1329       __os.fill(__space);
1330       __os.precision(std::numeric_limits<_RealType>::max_digits10);
1331 
1332       __os << __x.mu() << __space << __x.kappa();
1333 
1334       __os.flags(__flags);
1335       __os.fill(__fill);
1336       __os.precision(__precision);
1337       return __os;
1338     }
1339 
1340   template<typename _RealType, typename _CharT, typename _Traits>
1341     std::basic_istream<_CharT, _Traits>&
1342     operator>>(std::basic_istream<_CharT, _Traits>& __is,
1343 	       __gnu_cxx::von_mises_distribution<_RealType>& __x)
1344     {
1345       typedef std::basic_istream<_CharT, _Traits>  __istream_type;
1346       typedef typename __istream_type::ios_base    __ios_base;
1347 
1348       const typename __ios_base::fmtflags __flags = __is.flags();
1349       __is.flags(__ios_base::dec | __ios_base::skipws);
1350 
1351       _RealType __mu, __kappa;
1352       __is >> __mu >> __kappa;
1353       __x.param(typename __gnu_cxx::von_mises_distribution<_RealType>::
1354 		param_type(__mu, __kappa));
1355 
1356       __is.flags(__flags);
1357       return __is;
1358     }
1359 
1360 
1361   template<typename _UIntType>
1362     template<typename _UniformRandomNumberGenerator>
1363       typename hypergeometric_distribution<_UIntType>::result_type
1364       hypergeometric_distribution<_UIntType>::
1365       operator()(_UniformRandomNumberGenerator& __urng,
1366 		 const param_type& __param)
1367       {
1368 	std::__detail::_Adaptor<_UniformRandomNumberGenerator, double>
1369 	  __aurng(__urng);
1370 
1371 	result_type __a = __param.successful_size();
1372 	result_type __b = __param.total_size();
1373 	result_type __k = 0;
1374 
1375 	if (__param.total_draws() < __param.total_size() / 2)
1376 	  {
1377 	    for (result_type __i = 0; __i < __param.total_draws(); ++__i)
1378 	      {
1379 		if (__b * __aurng() < __a)
1380 		  {
1381 		    ++__k;
1382 		    if (__k == __param.successful_size())
1383 		      return __k;
1384 		   --__a;
1385 		  }
1386 		--__b;
1387 	      }
1388 	    return __k;
1389 	  }
1390 	else
1391 	  {
1392 	    for (result_type __i = 0; __i < __param.unsuccessful_size(); ++__i)
1393 	      {
1394 		if (__b * __aurng() < __a)
1395 		  {
1396 		    ++__k;
1397 		    if (__k == __param.successful_size())
1398 		      return __param.successful_size() - __k;
1399 		    --__a;
1400 		  }
1401 		--__b;
1402 	      }
1403 	    return __param.successful_size() - __k;
1404 	  }
1405       }
1406 
1407   template<typename _UIntType>
1408     template<typename _OutputIterator,
1409 	     typename _UniformRandomNumberGenerator>
1410       void
1411       hypergeometric_distribution<_UIntType>::
1412       __generate_impl(_OutputIterator __f, _OutputIterator __t,
1413 		      _UniformRandomNumberGenerator& __urng,
1414 		      const param_type& __param)
1415       {
1416 	__glibcxx_function_requires(_OutputIteratorConcept<_OutputIterator,
1417 	    result_type>)
1418 
1419 	while (__f != __t)
1420 	  *__f++ = this->operator()(__urng);
1421       }
1422 
1423   template<typename _UIntType, typename _CharT, typename _Traits>
1424     std::basic_ostream<_CharT, _Traits>&
1425     operator<<(std::basic_ostream<_CharT, _Traits>& __os,
1426 	       const __gnu_cxx::hypergeometric_distribution<_UIntType>& __x)
1427     {
1428       typedef std::basic_ostream<_CharT, _Traits>  __ostream_type;
1429       typedef typename __ostream_type::ios_base    __ios_base;
1430 
1431       const typename __ios_base::fmtflags __flags = __os.flags();
1432       const _CharT __fill = __os.fill();
1433       const std::streamsize __precision = __os.precision();
1434       const _CharT __space = __os.widen(' ');
1435       __os.flags(__ios_base::scientific | __ios_base::left);
1436       __os.fill(__space);
1437       __os.precision(std::numeric_limits<_UIntType>::max_digits10);
1438 
1439       __os << __x.total_size() << __space << __x.successful_size() << __space
1440 	   << __x.total_draws();
1441 
1442       __os.flags(__flags);
1443       __os.fill(__fill);
1444       __os.precision(__precision);
1445       return __os;
1446     }
1447 
1448   template<typename _UIntType, typename _CharT, typename _Traits>
1449     std::basic_istream<_CharT, _Traits>&
1450     operator>>(std::basic_istream<_CharT, _Traits>& __is,
1451 	       __gnu_cxx::hypergeometric_distribution<_UIntType>& __x)
1452     {
1453       typedef std::basic_istream<_CharT, _Traits>  __istream_type;
1454       typedef typename __istream_type::ios_base    __ios_base;
1455 
1456       const typename __ios_base::fmtflags __flags = __is.flags();
1457       __is.flags(__ios_base::dec | __ios_base::skipws);
1458 
1459       _UIntType __total_size, __successful_size, __total_draws;
1460       __is >> __total_size >> __successful_size >> __total_draws;
1461       __x.param(typename __gnu_cxx::hypergeometric_distribution<_UIntType>::
1462 		param_type(__total_size, __successful_size, __total_draws));
1463 
1464       __is.flags(__flags);
1465       return __is;
1466     }
1467 
1468 
1469   template<typename _RealType>
1470     template<typename _UniformRandomNumberGenerator>
1471       typename logistic_distribution<_RealType>::result_type
1472       logistic_distribution<_RealType>::
1473       operator()(_UniformRandomNumberGenerator& __urng,
1474 		 const param_type& __p)
1475       {
1476 	std::__detail::_Adaptor<_UniformRandomNumberGenerator, result_type>
1477 	  __aurng(__urng);
1478 
1479 	result_type __arg = result_type(1);
1480 	while (__arg == result_type(1) || __arg == result_type(0))
1481 	  __arg = __aurng();
1482 	return __p.a()
1483 	     + __p.b() * std::log(__arg / (result_type(1) - __arg));
1484       }
1485 
1486   template<typename _RealType>
1487     template<typename _OutputIterator,
1488 	     typename _UniformRandomNumberGenerator>
1489       void
1490       logistic_distribution<_RealType>::
1491       __generate_impl(_OutputIterator __f, _OutputIterator __t,
1492 		      _UniformRandomNumberGenerator& __urng,
1493 		      const param_type& __p)
1494       {
1495 	__glibcxx_function_requires(_OutputIteratorConcept<_OutputIterator,
1496 	    result_type>)
1497 
1498 	std::__detail::_Adaptor<_UniformRandomNumberGenerator, result_type>
1499 	  __aurng(__urng);
1500 
1501 	while (__f != __t)
1502 	  {
1503 	    result_type __arg = result_type(1);
1504 	    while (__arg == result_type(1) || __arg == result_type(0))
1505 	      __arg = __aurng();
1506 	    *__f++ = __p.a()
1507 		   + __p.b() * std::log(__arg / (result_type(1) - __arg));
1508 	  }
1509       }
1510 
1511   template<typename _RealType, typename _CharT, typename _Traits>
1512     std::basic_ostream<_CharT, _Traits>&
1513     operator<<(std::basic_ostream<_CharT, _Traits>& __os,
1514 	       const logistic_distribution<_RealType>& __x)
1515     {
1516       typedef std::basic_ostream<_CharT, _Traits>  __ostream_type;
1517       typedef typename __ostream_type::ios_base    __ios_base;
1518 
1519       const typename __ios_base::fmtflags __flags = __os.flags();
1520       const _CharT __fill = __os.fill();
1521       const std::streamsize __precision = __os.precision();
1522       const _CharT __space = __os.widen(' ');
1523       __os.flags(__ios_base::scientific | __ios_base::left);
1524       __os.fill(__space);
1525       __os.precision(std::numeric_limits<_RealType>::max_digits10);
1526 
1527       __os << __x.a() << __space << __x.b();
1528 
1529       __os.flags(__flags);
1530       __os.fill(__fill);
1531       __os.precision(__precision);
1532       return __os;
1533     }
1534 
1535   template<typename _RealType, typename _CharT, typename _Traits>
1536     std::basic_istream<_CharT, _Traits>&
1537     operator>>(std::basic_istream<_CharT, _Traits>& __is,
1538 	       logistic_distribution<_RealType>& __x)
1539     {
1540       typedef std::basic_istream<_CharT, _Traits>  __istream_type;
1541       typedef typename __istream_type::ios_base    __ios_base;
1542 
1543       const typename __ios_base::fmtflags __flags = __is.flags();
1544       __is.flags(__ios_base::dec | __ios_base::skipws);
1545 
1546       _RealType __a, __b;
1547       __is >> __a >> __b;
1548       __x.param(typename logistic_distribution<_RealType>::
1549 		param_type(__a, __b));
1550 
1551       __is.flags(__flags);
1552       return __is;
1553     }
1554 
1555 
1556   namespace {
1557 
1558     // Helper class for the uniform_on_sphere_distribution generation
1559     // function.
1560     template<std::size_t _Dimen, typename _RealType>
1561       class uniform_on_sphere_helper
1562       {
1563 	typedef typename uniform_on_sphere_distribution<_Dimen, _RealType>::
1564 	  result_type result_type;
1565 
1566       public:
1567 	template<typename _NormalDistribution,
1568 		 typename _UniformRandomNumberGenerator>
1569 	result_type operator()(_NormalDistribution& __nd,
1570 			       _UniformRandomNumberGenerator& __urng)
1571         {
1572 	  result_type __ret;
1573 	  typename result_type::value_type __norm;
1574 
1575 	  do
1576 	    {
1577 	      auto __sum = _RealType(0);
1578 
1579 	      std::generate(__ret.begin(), __ret.end(),
1580 			    [&__nd, &__urng, &__sum](){
1581 			      _RealType __t = __nd(__urng);
1582 			      __sum += __t * __t;
1583 			      return __t; });
1584 	      __norm = std::sqrt(__sum);
1585 	    }
1586 	  while (__norm == _RealType(0) || ! __builtin_isfinite(__norm));
1587 
1588 	  std::transform(__ret.begin(), __ret.end(), __ret.begin(),
1589 			 [__norm](_RealType __val){ return __val / __norm; });
1590 
1591 	  return __ret;
1592         }
1593       };
1594 
1595 
1596     template<typename _RealType>
1597       class uniform_on_sphere_helper<2, _RealType>
1598       {
1599 	typedef typename uniform_on_sphere_distribution<2, _RealType>::
1600 	  result_type result_type;
1601 
1602       public:
1603 	template<typename _NormalDistribution,
1604 		 typename _UniformRandomNumberGenerator>
1605 	result_type operator()(_NormalDistribution&,
1606 			       _UniformRandomNumberGenerator& __urng)
1607         {
1608 	  result_type __ret;
1609 	  _RealType __sq;
1610 	  std::__detail::_Adaptor<_UniformRandomNumberGenerator,
1611 				  _RealType> __aurng(__urng);
1612 
1613 	  do
1614 	    {
1615 	      __ret[0] = _RealType(2) * __aurng() - _RealType(1);
1616 	      __ret[1] = _RealType(2) * __aurng() - _RealType(1);
1617 
1618 	      __sq = __ret[0] * __ret[0] + __ret[1] * __ret[1];
1619 	    }
1620 	  while (__sq == _RealType(0) || __sq > _RealType(1));
1621 
1622 #if _GLIBCXX_USE_C99_MATH_TR1
1623 	  // Yes, we do not just use sqrt(__sq) because hypot() is more
1624 	  // accurate.
1625 	  auto __norm = std::hypot(__ret[0], __ret[1]);
1626 #else
1627 	  auto __norm = std::sqrt(__sq);
1628 #endif
1629 	  __ret[0] /= __norm;
1630 	  __ret[1] /= __norm;
1631 
1632 	  return __ret;
1633         }
1634       };
1635 
1636   }
1637 
1638 
1639   template<std::size_t _Dimen, typename _RealType>
1640     template<typename _UniformRandomNumberGenerator>
1641       typename uniform_on_sphere_distribution<_Dimen, _RealType>::result_type
1642       uniform_on_sphere_distribution<_Dimen, _RealType>::
1643       operator()(_UniformRandomNumberGenerator& __urng,
1644 		 const param_type& __p)
1645       {
1646         uniform_on_sphere_helper<_Dimen, _RealType> __helper;
1647         return __helper(_M_nd, __urng);
1648       }
1649 
1650   template<std::size_t _Dimen, typename _RealType>
1651     template<typename _OutputIterator,
1652 	     typename _UniformRandomNumberGenerator>
1653       void
1654       uniform_on_sphere_distribution<_Dimen, _RealType>::
1655       __generate_impl(_OutputIterator __f, _OutputIterator __t,
1656 		      _UniformRandomNumberGenerator& __urng,
1657 		      const param_type& __param)
1658       {
1659 	__glibcxx_function_requires(_OutputIteratorConcept<_OutputIterator,
1660 	    result_type>)
1661 
1662 	while (__f != __t)
1663 	  *__f++ = this->operator()(__urng, __param);
1664       }
1665 
1666   template<std::size_t _Dimen, typename _RealType, typename _CharT,
1667 	   typename _Traits>
1668     std::basic_ostream<_CharT, _Traits>&
1669     operator<<(std::basic_ostream<_CharT, _Traits>& __os,
1670 	       const __gnu_cxx::uniform_on_sphere_distribution<_Dimen,
1671 							       _RealType>& __x)
1672     {
1673       return __os << __x._M_nd;
1674     }
1675 
1676   template<std::size_t _Dimen, typename _RealType, typename _CharT,
1677 	   typename _Traits>
1678     std::basic_istream<_CharT, _Traits>&
1679     operator>>(std::basic_istream<_CharT, _Traits>& __is,
1680 	       __gnu_cxx::uniform_on_sphere_distribution<_Dimen,
1681 							 _RealType>& __x)
1682     {
1683       return __is >> __x._M_nd;
1684     }
1685 
1686 
1687   namespace {
1688 
1689     // Helper class for the uniform_inside_sphere_distribution generation
1690     // function.
1691     template<std::size_t _Dimen, bool _SmallDimen, typename _RealType>
1692       class uniform_inside_sphere_helper;
1693 
1694     template<std::size_t _Dimen, typename _RealType>
1695       class uniform_inside_sphere_helper<_Dimen, false, _RealType>
1696       {
1697 	using result_type
1698 	  = typename uniform_inside_sphere_distribution<_Dimen, _RealType>::
1699 	    result_type;
1700 
1701       public:
1702 	template<typename _UniformOnSphereDistribution,
1703 		 typename _UniformRandomNumberGenerator>
1704 	result_type
1705 	operator()(_UniformOnSphereDistribution& __uosd,
1706 		   _UniformRandomNumberGenerator& __urng,
1707 		   _RealType __radius)
1708         {
1709 	  std::__detail::_Adaptor<_UniformRandomNumberGenerator,
1710 				  _RealType> __aurng(__urng);
1711 
1712 	  _RealType __pow = 1 / _RealType(_Dimen);
1713 	  _RealType __urt = __radius * std::pow(__aurng(), __pow);
1714 	  result_type __ret = __uosd(__aurng);
1715 
1716 	  std::transform(__ret.begin(), __ret.end(), __ret.begin(),
1717 			 [__urt](_RealType __val)
1718 			 { return __val * __urt; });
1719 
1720 	  return __ret;
1721         }
1722       };
1723 
1724     // Helper class for the uniform_inside_sphere_distribution generation
1725     // function specialized for small dimensions.
1726     template<std::size_t _Dimen, typename _RealType>
1727       class uniform_inside_sphere_helper<_Dimen, true, _RealType>
1728       {
1729 	using result_type
1730 	  = typename uniform_inside_sphere_distribution<_Dimen, _RealType>::
1731 	    result_type;
1732 
1733       public:
1734 	template<typename _UniformOnSphereDistribution,
1735 		 typename _UniformRandomNumberGenerator>
1736 	result_type
1737 	operator()(_UniformOnSphereDistribution&,
1738 		   _UniformRandomNumberGenerator& __urng,
1739 		   _RealType __radius)
1740         {
1741 	  result_type __ret;
1742 	  _RealType __sq;
1743 	  _RealType __radsq = __radius * __radius;
1744 	  std::__detail::_Adaptor<_UniformRandomNumberGenerator,
1745 				  _RealType> __aurng(__urng);
1746 
1747 	  do
1748 	    {
1749 	      __sq = _RealType(0);
1750 	      for (int i = 0; i < _Dimen; ++i)
1751 		{
1752 		  __ret[i] = _RealType(2) * __aurng() - _RealType(1);
1753 		  __sq += __ret[i] * __ret[i];
1754 		}
1755 	    }
1756 	  while (__sq > _RealType(1));
1757 
1758 	  for (int i = 0; i < _Dimen; ++i)
1759             __ret[i] *= __radius;
1760 
1761 	  return __ret;
1762         }
1763       };
1764   } // namespace
1765 
1766   //
1767   //  Experiments have shown that rejection is more efficient than transform
1768   //  for dimensions less than 8.
1769   //
1770   template<std::size_t _Dimen, typename _RealType>
1771     template<typename _UniformRandomNumberGenerator>
1772       typename uniform_inside_sphere_distribution<_Dimen, _RealType>::result_type
1773       uniform_inside_sphere_distribution<_Dimen, _RealType>::
1774       operator()(_UniformRandomNumberGenerator& __urng,
1775 		 const param_type& __p)
1776       {
1777         uniform_inside_sphere_helper<_Dimen, _Dimen < 8, _RealType> __helper;
1778         return __helper(_M_uosd, __urng, __p.radius());
1779       }
1780 
1781   template<std::size_t _Dimen, typename _RealType>
1782     template<typename _OutputIterator,
1783 	     typename _UniformRandomNumberGenerator>
1784       void
1785       uniform_inside_sphere_distribution<_Dimen, _RealType>::
1786       __generate_impl(_OutputIterator __f, _OutputIterator __t,
1787 		      _UniformRandomNumberGenerator& __urng,
1788 		      const param_type& __param)
1789       {
1790 	__glibcxx_function_requires(_OutputIteratorConcept<_OutputIterator,
1791 	    result_type>)
1792 
1793 	while (__f != __t)
1794 	  *__f++ = this->operator()(__urng, __param);
1795       }
1796 
1797   template<std::size_t _Dimen, typename _RealType, typename _CharT,
1798 	   typename _Traits>
1799     std::basic_ostream<_CharT, _Traits>&
1800     operator<<(std::basic_ostream<_CharT, _Traits>& __os,
1801 	       const __gnu_cxx::uniform_inside_sphere_distribution<_Dimen,
1802 								_RealType>& __x)
1803     {
1804       typedef std::basic_ostream<_CharT, _Traits>  __ostream_type;
1805       typedef typename __ostream_type::ios_base    __ios_base;
1806 
1807       const typename __ios_base::fmtflags __flags = __os.flags();
1808       const _CharT __fill = __os.fill();
1809       const std::streamsize __precision = __os.precision();
1810       const _CharT __space = __os.widen(' ');
1811       __os.flags(__ios_base::scientific | __ios_base::left);
1812       __os.fill(__space);
1813       __os.precision(std::numeric_limits<_RealType>::max_digits10);
1814 
1815       __os << __x.radius() << __space << __x._M_uosd;
1816 
1817       __os.flags(__flags);
1818       __os.fill(__fill);
1819       __os.precision(__precision);
1820 
1821       return __os;
1822     }
1823 
1824   template<std::size_t _Dimen, typename _RealType, typename _CharT,
1825 	   typename _Traits>
1826     std::basic_istream<_CharT, _Traits>&
1827     operator>>(std::basic_istream<_CharT, _Traits>& __is,
1828 	       __gnu_cxx::uniform_inside_sphere_distribution<_Dimen,
1829 							     _RealType>& __x)
1830     {
1831       typedef std::basic_istream<_CharT, _Traits>  __istream_type;
1832       typedef typename __istream_type::ios_base    __ios_base;
1833 
1834       const typename __ios_base::fmtflags __flags = __is.flags();
1835       __is.flags(__ios_base::dec | __ios_base::skipws);
1836 
1837       _RealType __radius_val;
1838       __is >> __radius_val >> __x._M_uosd;
1839       __x.param(typename uniform_inside_sphere_distribution<_Dimen, _RealType>::
1840 		param_type(__radius_val));
1841 
1842       __is.flags(__flags);
1843 
1844       return __is;
1845     }
1846 
1847 _GLIBCXX_END_NAMESPACE_VERSION
1848 } // namespace __gnu_cxx
1849 
1850 
1851 #endif // _EXT_RANDOM_TCC
1852