blob: dac6dd11703581e2a3a3b16fa72fdebacea86cb4 [file] [log] [blame]
Howard Hinnantbc8d3f92010-05-11 19:42:16 +00001// -*- C++ -*-
2//===--------------------------- random -----------------------------------===//
3//
Howard Hinnantf5256e12010-05-11 21:36:01 +00004// The LLVM Compiler Infrastructure
Howard Hinnantbc8d3f92010-05-11 19:42:16 +00005//
6// This file is distributed under the University of Illinois Open Source
7// License. See LICENSE.TXT for details.
8//
9//===----------------------------------------------------------------------===//
10
11#ifndef _LIBCPP_RANDOM
12#define _LIBCPP_RANDOM
13
14/*
15 random synopsis
16
17#include <initializer_list>
18
19namespace std
20{
21
22// Engines
23
24template <class UIntType, UIntType a, UIntType c, UIntType m>
25class linear_congruential_engine
26{
27public:
28 // types
29 typedef UIntType result_type;
30
31 // engine characteristics
32 static constexpr result_type multiplier = a;
33 static constexpr result_type increment = c;
34 static constexpr result_type modulus = m;
35 static constexpr result_type min() { return c == 0u ? 1u: 0u;}
36 static constexpr result_type max() { return m - 1u;}
37 static constexpr result_type default_seed = 1u;
38
39 // constructors and seeding functions
40 explicit linear_congruential_engine(result_type s = default_seed);
41 template<class Sseq> explicit linear_congruential_engine(Sseq& q);
42 void seed(result_type s = default_seed);
43 template<class Sseq> void seed(Sseq& q);
44
45 // generating functions
46 result_type operator()();
47 void discard(unsigned long long z);
48};
49
50template <class UIntType, UIntType a, UIntType c, UIntType m>
51bool
52operator==(const linear_congruential_engine<UIntType, a, c, m>& x,
53 const linear_congruential_engine<UIntType, a, c, m>& y);
54
55template <class UIntType, UIntType a, UIntType c, UIntType m>
56bool
57operator!=(const linear_congruential_engine<UIntType, a, c, m>& x,
58 const linear_congruential_engine<UIntType, a, c, m>& y);
59
60template <class charT, class traits,
61 class UIntType, UIntType a, UIntType c, UIntType m>
62basic_ostream<charT, traits>&
63operator<<(basic_ostream<charT, traits>& os,
64 const linear_congruential_engine<UIntType, a, c, m>& x);
65
66template <class charT, class traits,
67 class UIntType, UIntType a, UIntType c, UIntType m>
68basic_istream<charT, traits>&
69operator>>(basic_istream<charT, traits>& is,
70 linear_congruential_engine<UIntType, a, c, m>& x);
71
72template <class UIntType, size_t w, size_t n, size_t m, size_t r,
73 UIntType a, size_t u, UIntType d, size_t s,
74 UIntType b, size_t t, UIntType c, size_t l, UIntType f>
75class mersenne_twister_engine
76{
77public:
78 // types
79 typedef UIntType result_type;
80
81 // engine characteristics
82 static constexpr size_t word_size = w;
83 static constexpr size_t state_size = n;
84 static constexpr size_t shift_size = m;
85 static constexpr size_t mask_bits = r;
86 static constexpr result_type xor_mask = a;
87 static constexpr size_t tempering_u = u;
88 static constexpr result_type tempering_d = d;
89 static constexpr size_t tempering_s = s;
90 static constexpr result_type tempering_b = b;
91 static constexpr size_t tempering_t = t;
92 static constexpr result_type tempering_c = c;
93 static constexpr size_t tempering_l = l;
94 static constexpr result_type initialization_multiplier = f;
95 static constexpr result_type min () { return 0; }
96 static constexpr result_type max() { return 2^w - 1; }
97 static constexpr result_type default_seed = 5489u;
98
99 // constructors and seeding functions
100 explicit mersenne_twister_engine(result_type value = default_seed);
101 template<class Sseq> explicit mersenne_twister_engine(Sseq& q);
102 void seed(result_type value = default_seed);
103 template<class Sseq> void seed(Sseq& q);
104
105 // generating functions
106 result_type operator()();
107 void discard(unsigned long long z);
108};
109
110template <class UIntType, size_t w, size_t n, size_t m, size_t r,
111 UIntType a, size_t u, UIntType d, size_t s,
112 UIntType b, size_t t, UIntType c, size_t l, UIntType f>
113bool
114operator==(
115 const mersenne_twister_engine<UIntType, w, n, m, r, a, u, d, s, b, t, c, l, f>& x,
116 const mersenne_twister_engine<UIntType, w, n, m, r, a, u, d, s, b, t, c, l, f>& y);
117
118template <class UIntType, size_t w, size_t n, size_t m, size_t r,
119 UIntType a, size_t u, UIntType d, size_t s,
120 UIntType b, size_t t, UIntType c, size_t l, UIntType f>
121bool
122operator!=(
123 const mersenne_twister_engine<UIntType, w, n, m, r, a, u, d, s, b, t, c, l, f>& x,
124 const mersenne_twister_engine<UIntType, w, n, m, r, a, u, d, s, b, t, c, l, f>& y);
125
126template <class charT, class traits,
127 class UIntType, size_t w, size_t n, size_t m, size_t r,
128 UIntType a, size_t u, UIntType d, size_t s,
129 UIntType b, size_t t, UIntType c, size_t l, UIntType f>
130basic_ostream<charT, traits>&
131operator<<(basic_ostream<charT, traits>& os,
132 const mersenne_twister_engine<UIntType, w, n, m, r, a, u, d, s, b, t, c, l, f>& x);
133
134template <class charT, class traits,
135 class UIntType, size_t w, size_t n, size_t m, size_t r,
136 UIntType a, size_t u, UIntType d, size_t s,
137 UIntType b, size_t t, UIntType c, size_t l, UIntType f>
138basic_istream<charT, traits>&
139operator>>(basic_istream<charT, traits>& is,
140 mersenne_twister_engine<UIntType, w, n, m, r, a, u, d, s, b, t, c, l, f>& x);
141
142template<class UIntType, size_t w, size_t s, size_t r>
143class subtract_with_carry_engine
144{
145public:
146 // types
147 typedef UIntType result_type;
148
149 // engine characteristics
150 static constexpr size_t word_size = w;
151 static constexpr size_t short_lag = s;
152 static constexpr size_t long_lag = r;
153 static constexpr result_type min() { return 0; }
154 static constexpr result_type max() { return m-1; }
155 static constexpr result_type default_seed = 19780503u;
156
157 // constructors and seeding functions
158 explicit subtract_with_carry_engine(result_type value = default_seed);
159 template<class Sseq> explicit subtract_with_carry_engine(Sseq& q);
160 void seed(result_type value = default_seed);
161 template<class Sseq> void seed(Sseq& q);
162
163 // generating functions
164 result_type operator()();
165 void discard(unsigned long long z);
166};
167
168template<class UIntType, size_t w, size_t s, size_t r>
169bool
170operator==(
171 const subtract_with_carry_engine<UIntType, w, s, r>& x,
172 const subtract_with_carry_engine<UIntType, w, s, r>& y);
173
174template<class UIntType, size_t w, size_t s, size_t r>
175bool
176operator!=(
177 const subtract_with_carry_engine<UIntType, w, s, r>& x,
178 const subtract_with_carry_engine<UIntType, w, s, r>& y);
179
180template <class charT, class traits,
181 class UIntType, size_t w, size_t s, size_t r>
182basic_ostream<charT, traits>&
183operator<<(basic_ostream<charT, traits>& os,
184 const subtract_with_carry_engine<UIntType, w, s, r>& x);
185
186template <class charT, class traits,
187 class UIntType, size_t w, size_t s, size_t r>
188basic_istream<charT, traits>&
189operator>>(basic_istream<charT, traits>& is,
190 subtract_with_carry_engine<UIntType, w, s, r>& x);
191
192template<class Engine, size_t p, size_t r>
193class discard_block_engine
194{
195public:
196 // types
197 typedef typename Engine::result_type result_type;
198
199 // engine characteristics
200 static constexpr size_t block_size = p;
201 static constexpr size_t used_block = r;
202 static constexpr result_type min() { return Engine::min(); }
203 static constexpr result_type max() { return Engine::max(); }
204
205 // constructors and seeding functions
206 discard_block_engine();
207 explicit discard_block_engine(const Engine& e);
208 explicit discard_block_engine(Engine&& e);
209 explicit discard_block_engine(result_type s);
210 template<class Sseq> explicit discard_block_engine(Sseq& q);
211 void seed();
212 void seed(result_type s);
213 template<class Sseq> void seed(Sseq& q);
214
215 // generating functions
216 result_type operator()();
217 void discard(unsigned long long z);
218
219 // property functions
220 const Engine& base() const;
221};
222
223template<class Engine, size_t p, size_t r>
224bool
225operator==(
226 const discard_block_engine<Engine, p, r>& x,
227 const discard_block_engine<Engine, p, r>& y);
228
229template<class Engine, size_t p, size_t r>
230bool
231operator!=(
232 const discard_block_engine<Engine, p, r>& x,
233 const discard_block_engine<Engine, p, r>& y);
234
235template <class charT, class traits,
236 class Engine, size_t p, size_t r>
237basic_ostream<charT, traits>&
238operator<<(basic_ostream<charT, traits>& os,
239 const discard_block_engine<Engine, p, r>& x);
240
241template <class charT, class traits,
242 class Engine, size_t p, size_t r>
243basic_istream<charT, traits>&
244operator>>(basic_istream<charT, traits>& is,
245 discard_block_engine<Engine, p, r>& x);
246
247template<class Engine, size_t w, class UIntType>
248class independent_bits_engine
249{
250public:
251 // types
252 typedef UIntType result_type;
253
254 // engine characteristics
255 static constexpr result_type min() { return 0; }
256 static constexpr result_type max() { return 2^w - 1; }
257
258 // constructors and seeding functions
259 independent_bits_engine();
260 explicit independent_bits_engine(const Engine& e);
261 explicit independent_bits_engine(Engine&& e);
262 explicit independent_bits_engine(result_type s);
263 template<class Sseq> explicit independent_bits_engine(Sseq& q);
264 void seed();
265 void seed(result_type s);
266 template<class Sseq> void seed(Sseq& q);
267
268 // generating functions
269 result_type operator()(); void discard(unsigned long long z);
270
271 // property functions
272 const Engine& base() const;
273};
274
275template<class Engine, size_t w, class UIntType>
276bool
277operator==(
278 const independent_bits_engine<Engine, w, UIntType>& x,
279 const independent_bits_engine<Engine, w, UIntType>& y);
280
281template<class Engine, size_t w, class UIntType>
282bool
283operator!=(
284 const independent_bits_engine<Engine, w, UIntType>& x,
285 const independent_bits_engine<Engine, w, UIntType>& y);
286
287template <class charT, class traits,
288 class Engine, size_t w, class UIntType>
289basic_ostream<charT, traits>&
290operator<<(basic_ostream<charT, traits>& os,
291 const independent_bits_engine<Engine, w, UIntType>& x);
292
293template <class charT, class traits,
294 class Engine, size_t w, class UIntType>
295basic_istream<charT, traits>&
296operator>>(basic_istream<charT, traits>& is,
297 independent_bits_engine<Engine, w, UIntType>& x);
298
299template<class Engine, size_t k>
300class shuffle_order_engine
301{
302public:
303 // types
304 typedef typename Engine::result_type result_type;
305
306 // engine characteristics
307 static constexpr size_t table_size = k;
308 static constexpr result_type min() { return Engine::min; }
309 static constexpr result_type max() { return Engine::max; }
310
311 // constructors and seeding functions
312 shuffle_order_engine();
313 explicit shuffle_order_engine(const Engine& e);
314 explicit shuffle_order_engine(Engine&& e);
315 explicit shuffle_order_engine(result_type s);
316 template<class Sseq> explicit shuffle_order_engine(Sseq& q);
317 void seed();
318 void seed(result_type s);
319 template<class Sseq> void seed(Sseq& q);
320
321 // generating functions
322 result_type operator()();
323 void discard(unsigned long long z);
324
325 // property functions
326 const Engine& base() const;
327};
328
329template<class Engine, size_t k>
330bool
331operator==(
332 const shuffle_order_engine<Engine, k>& x,
333 const shuffle_order_engine<Engine, k>& y);
334
335template<class Engine, size_t k>
336bool
337operator!=(
338 const shuffle_order_engine<Engine, k>& x,
339 const shuffle_order_engine<Engine, k>& y);
340
341template <class charT, class traits,
342 class Engine, size_t k>
343basic_ostream<charT, traits>&
344operator<<(basic_ostream<charT, traits>& os,
345 const shuffle_order_engine<Engine, k>& x);
346
347template <class charT, class traits,
348 class Engine, size_t k>
349basic_istream<charT, traits>&
350operator>>(basic_istream<charT, traits>& is,
351 shuffle_order_engine<Engine, k>& x);
352
353typedef linear_congruential_engine<uint_fast32_t, 16807, 0, 2147483647>
354 minstd_rand0;
355typedef linear_congruential_engine<uint_fast32_t, 48271, 0, 2147483647>
356 minstd_rand;
357typedef mersenne_twister_engine<uint_fast32_t, 32, 624, 397, 31,
358 0x9908b0df,
359 11, 0xffffffff,
360 7, 0x9d2c5680,
361 15, 0xefc60000,
362 18, 1812433253> mt19937;
363typedef mersenne_twister_engine<uint_fast64_t, 64, 312, 156, 31,
364 0xb5026f5aa96619e9,
365 29, 0x5555555555555555,
366 17, 0x71d67fffeda60000,
367 37, 0xfff7eee000000000,
368 43, 6364136223846793005> mt19937_64;
369typedef subtract_with_carry_engine<uint_fast32_t, 24, 10, 24> ranlux24_base;
370typedef subtract_with_carry_engine<uint_fast64_t, 48, 5, 12> ranlux48_base;
371typedef discard_block_engine<ranlux24_base, 223, 23> ranlux24;
372typedef discard_block_engine<ranlux48_base, 389, 11> ranlux48;
373typedef shuffle_order_engine<minstd_rand0, 256> knuth_b;
374typedef minstd_rand0 default_random_engine;
375
376// Generators
377
378class random_device
379{
380public:
381 // types
382 typedef unsigned int result_type;
383
384 // generator characteristics
385 static constexpr result_type min() { return numeric_limits<result_type>::min(); }
386 static constexpr result_type max() { return numeric_limits<result_type>::max(); }
387
388 // constructors
389 explicit random_device(const string& token = "/dev/urandom");
390
391 // generating functions
392 result_type operator()();
393
394 // property functions
395 double entropy() const;
396
397 // no copy functions
398 random_device(const random_device& ) = delete;
399 void operator=(const random_device& ) = delete;
400};
401
402// Utilities
403
404class seed_seq
405{
406public:
407 // types
408 typedef uint_least32_t result_type;
409
410 // constructors
411 seed_seq();
412 template<class T>
413 seed_seq(initializer_list<T> il);
414 template<class InputIterator>
415 seed_seq(InputIterator begin, InputIterator end);
416
417 // generating functions
418 template<class RandomAccessIterator>
419 void generate(RandomAccessIterator begin, RandomAccessIterator end);
420
421 // property functions
422 size_t size() const;
423 template<class OutputIterator>
424 void param(OutputIterator dest) const;
425
426 // no copy functions
427 seed_seq(const seed_seq&) = delete;
428 void operator=(const seed_seq& ) = delete;
429};
430
431template<class RealType, size_t bits, class URNG>
432 RealType generate_canonical(URNG& g);
433
434// Distributions
435
436template<class IntType = int>
437class uniform_int_distribution
438{
439public:
440 // types
441 typedef IntType result_type;
442
443 class param_type
444 {
445 public:
446 typedef uniform_int_distribution distribution_type;
447
448 explicit param_type(IntType a = 0,
449 IntType b = numeric_limits<IntType>::max());
450
451 result_type a() const;
452 result_type b() const;
453
454 friend bool operator==(const param_type& x, const param_type& y);
455 friend bool operator!=(const param_type& x, const param_type& y);
456 };
457
458 // constructors and reset functions
459 explicit uniform_int_distribution(IntType a = 0,
460 IntType b = numeric_limits<IntType>::max());
461 explicit uniform_int_distribution(const param_type& parm);
462 void reset();
463
464 // generating functions
465 template<class URNG> result_type operator()(URNG& g);
466 template<class URNG> result_type operator()(URNG& g, const param_type& parm);
467
468 // property functions
469 result_type a() const;
470 result_type b() const;
471
472 param_type param() const;
473 void param(const param_type& parm);
474
475 result_type min() const;
476 result_type max() const;
477
478 friend bool operator==(const uniform_int_distribution& x,
479 const uniform_int_distribution& y);
480 friend bool operator!=(const uniform_int_distribution& x,
481 const uniform_int_distribution& y);
482
483 template <class charT, class traits>
484 friend
485 basic_ostream<charT, traits>&
486 operator<<(basic_ostream<charT, traits>& os,
487 const uniform_int_distribution& x);
488
489 template <class charT, class traits>
490 friend
491 basic_istream<charT, traits>&
492 operator>>(basic_istream<charT, traits>& is,
493 uniform_int_distribution& x);
494};
495
496template<class RealType = double>
497class uniform_real_distribution
498{
499public:
500 // types
501 typedef RealType result_type;
502
503 class param_type
504 {
505 public:
506 typedef uniform_real_distribution distribution_type;
507
508 explicit param_type(RealType a = 0,
509 RealType b = 1);
510
511 result_type a() const;
512 result_type b() const;
513
514 friend bool operator==(const param_type& x, const param_type& y);
515 friend bool operator!=(const param_type& x, const param_type& y);
516 };
517
518 // constructors and reset functions
519 explicit uniform_real_distribution(RealType a = 0.0, RealType b = 1.0);
520 explicit uniform_real_distribution(const param_type& parm);
521 void reset();
522
523 // generating functions
524 template<class URNG> result_type operator()(URNG& g);
525 template<class URNG> result_type operator()(URNG& g, const param_type& parm);
526
527 // property functions
528 result_type a() const;
529 result_type b() const;
530
531 param_type param() const;
532 void param(const param_type& parm);
533
534 result_type min() const;
535 result_type max() const;
536
537 friend bool operator==(const uniform_real_distribution& x,
538 const uniform_real_distribution& y);
539 friend bool operator!=(const uniform_real_distribution& x,
540 const uniform_real_distribution& y);
541
542 template <class charT, class traits>
543 friend
544 basic_ostream<charT, traits>&
545 operator<<(basic_ostream<charT, traits>& os,
546 const uniform_real_distribution& x);
547
548 template <class charT, class traits>
549 friend
550 basic_istream<charT, traits>&
551 operator>>(basic_istream<charT, traits>& is,
552 uniform_real_distribution& x);
553};
554
555class bernoulli_distribution
556{
557public:
558 // types
559 typedef bool result_type;
560
561 class param_type
562 {
563 public:
564 typedef bernoulli_distribution distribution_type;
565
566 explicit param_type(double p = 0.5);
567
568 double p() const;
569
570 friend bool operator==(const param_type& x, const param_type& y);
571 friend bool operator!=(const param_type& x, const param_type& y);
572 };
573
574 // constructors and reset functions
575 explicit bernoulli_distribution(double p = 0.5);
576 explicit bernoulli_distribution(const param_type& parm);
577 void reset();
578
579 // generating functions
580 template<class URNG> result_type operator()(URNG& g);
581 template<class URNG> result_type operator()(URNG& g, const param_type& parm);
582
583 // property functions
584 double p() const;
585
586 param_type param() const;
587 void param(const param_type& parm);
588
589 result_type min() const;
590 result_type max() const;
591
592 friend bool operator==(const bernoulli_distribution& x,
593 const bernoulli_distribution& y);
594 friend bool operator!=(const bernoulli_distribution& x,
595 const bernoulli_distribution& y);
596
597 template <class charT, class traits>
598 friend
599 basic_ostream<charT, traits>&
600 operator<<(basic_ostream<charT, traits>& os,
601 const bernoulli_distribution& x);
602
603 template <class charT, class traits>
604 friend
605 basic_istream<charT, traits>&
606 operator>>(basic_istream<charT, traits>& is,
607 bernoulli_distribution& x);
608};
609
610template<class IntType = int>
Howard Hinnant03aad812010-05-11 23:26:59 +0000611class binomial_distribution
612{
613public:
614 // types
615 typedef IntType result_type;
616
617 class param_type
618 {
619 public:
620 typedef binomial_distribution distribution_type;
621
622 explicit param_type(IntType t = 1, double p = 0.5);
623
624 IntType t() const;
625 double p() const;
626
627 friend bool operator==(const param_type& x, const param_type& y);
628 friend bool operator!=(const param_type& x, const param_type& y);
629 };
630
631 // constructors and reset functions
632 explicit binomial_distribution(IntType t = 1, double p = 0.5);
633 explicit binomial_distribution(const param_type& parm);
634 void reset();
635
636 // generating functions
637 template<class URNG> result_type operator()(URNG& g);
638 template<class URNG> result_type operator()(URNG& g, const param_type& parm);
639
640 // property functions
641 IntType t() const;
642 double p() const;
643
644 param_type param() const;
645 void param(const param_type& parm);
646
647 result_type min() const;
648 result_type max() const;
649
650 friend bool operator==(const binomial_distribution& x,
651 const binomial_distribution& y);
652 friend bool operator!=(const binomial_distribution& x,
653 const binomial_distribution& y);
654
655 template <class charT, class traits>
656 friend
657 basic_ostream<charT, traits>&
658 operator<<(basic_ostream<charT, traits>& os,
659 const binomial_distribution& x);
660
661 template <class charT, class traits>
662 friend
663 basic_istream<charT, traits>&
664 operator>>(basic_istream<charT, traits>& is,
665 binomial_distribution& x);
666};
Howard Hinnantbc8d3f92010-05-11 19:42:16 +0000667
668template<class IntType = int>
669 class geometric_distribution;
670
671template<class IntType = int>
672 class negative_binomial_distribution;
673
674template<class IntType = int>
Howard Hinnant4ff556c2010-05-14 21:38:54 +0000675class poisson_distribution
676{
677public:
678 // types
679 typedef IntType result_type;
680
681 class param_type
682 {
683 public:
684 typedef poisson_distribution distribution_type;
685
686 explicit param_type(double mean = 1.0);
687
688 double mean() const;
689
690 friend bool operator==(const param_type& x, const param_type& y);
691 friend bool operator!=(const param_type& x, const param_type& y);
692 };
693
694 // constructors and reset functions
695 explicit poisson_distribution(double mean = 1.0);
696 explicit poisson_distribution(const param_type& parm);
697 void reset();
698
699 // generating functions
700 template<class URNG> result_type operator()(URNG& g);
701 template<class URNG> result_type operator()(URNG& g, const param_type& parm);
702
703 // property functions
704 double mean() const;
705
706 param_type param() const;
707 void param(const param_type& parm);
708
709 result_type min() const;
710 result_type max() const;
711
712 friend bool operator==(const poisson_distribution& x,
713 const poisson_distribution& y);
714 friend bool operator!=(const poisson_distribution& x,
715 const poisson_distribution& y);
716
717 template <class charT, class traits>
718 friend
719 basic_ostream<charT, traits>&
720 operator<<(basic_ostream<charT, traits>& os,
721 const poisson_distribution& x);
722
723 template <class charT, class traits>
724 friend
725 basic_istream<charT, traits>&
726 operator>>(basic_istream<charT, traits>& is,
727 poisson_distribution& x);
728};
Howard Hinnantbc8d3f92010-05-11 19:42:16 +0000729
730template<class RealType = double>
Howard Hinnant30a840f2010-05-12 17:08:57 +0000731class exponential_distribution
732{
733public:
734 // types
735 typedef RealType result_type;
736
737 class param_type
738 {
739 public:
740 typedef exponential_distribution distribution_type;
741
Howard Hinnanta64111c2010-05-12 21:02:31 +0000742 explicit param_type(result_type lambda = 1.0);
Howard Hinnant30a840f2010-05-12 17:08:57 +0000743
Howard Hinnanta64111c2010-05-12 21:02:31 +0000744 result_type lambda() const;
Howard Hinnant30a840f2010-05-12 17:08:57 +0000745
746 friend bool operator==(const param_type& x, const param_type& y);
747 friend bool operator!=(const param_type& x, const param_type& y);
748 };
749
750 // constructors and reset functions
Howard Hinnanta64111c2010-05-12 21:02:31 +0000751 explicit exponential_distribution(result_type lambda = 1.0);
Howard Hinnant30a840f2010-05-12 17:08:57 +0000752 explicit exponential_distribution(const param_type& parm);
753 void reset();
754
755 // generating functions
756 template<class URNG> result_type operator()(URNG& g);
757 template<class URNG> result_type operator()(URNG& g, const param_type& parm);
758
759 // property functions
Howard Hinnanta64111c2010-05-12 21:02:31 +0000760 result_type lambda() const;
Howard Hinnant30a840f2010-05-12 17:08:57 +0000761
762 param_type param() const;
763 void param(const param_type& parm);
764
765 result_type min() const;
766 result_type max() const;
767
768 friend bool operator==(const exponential_distribution& x,
769 const exponential_distribution& y);
770 friend bool operator!=(const exponential_distribution& x,
771 const exponential_distribution& y);
772
773 template <class charT, class traits>
774 friend
775 basic_ostream<charT, traits>&
776 operator<<(basic_ostream<charT, traits>& os,
777 const exponential_distribution& x);
778
779 template <class charT, class traits>
780 friend
781 basic_istream<charT, traits>&
782 operator>>(basic_istream<charT, traits>& is,
783 exponential_distribution& x);
784};
Howard Hinnantbc8d3f92010-05-11 19:42:16 +0000785
786template<class RealType = double>
Howard Hinnantc7c49132010-05-13 17:58:28 +0000787class gamma_distribution
788{
789public:
790 // types
791 typedef RealType result_type;
792
793 class param_type
794 {
795 public:
796 typedef gamma_distribution distribution_type;
797
798 explicit param_type(result_type alpha = 1, result_type beta = 1);
799
800 result_type alpha() const;
801 result_type beta() const;
802
803 friend bool operator==(const param_type& x, const param_type& y);
804 friend bool operator!=(const param_type& x, const param_type& y);
805 };
806
807 // constructors and reset functions
808 explicit gamma_distribution(result_type alpha = 1, result_type beta = 1);
809 explicit gamma_distribution(const param_type& parm);
810 void reset();
811
812 // generating functions
813 template<class URNG> result_type operator()(URNG& g);
814 template<class URNG> result_type operator()(URNG& g, const param_type& parm);
815
816 // property functions
817 result_type alpha() const;
818 result_type beta() const;
819
820 param_type param() const;
821 void param(const param_type& parm);
822
823 result_type min() const;
824 result_type max() const;
825
826 friend bool operator==(const gamma_distribution& x,
827 const gamma_distribution& y);
828 friend bool operator!=(const gamma_distribution& x,
829 const gamma_distribution& y);
830
831 template <class charT, class traits>
832 friend
833 basic_ostream<charT, traits>&
834 operator<<(basic_ostream<charT, traits>& os,
835 const gamma_distribution& x);
836
837 template <class charT, class traits>
838 friend
839 basic_istream<charT, traits>&
840 operator>>(basic_istream<charT, traits>& is,
841 gamma_distribution& x);
842};
Howard Hinnantbc8d3f92010-05-11 19:42:16 +0000843
844template<class RealType = double>
Howard Hinnant9de6e302010-05-16 01:09:02 +0000845class weibull_distribution
846{
847public:
848 // types
849 typedef RealType result_type;
850
851 class param_type
852 {
853 public:
854 typedef weibull_distribution distribution_type;
855
856 explicit param_type(result_type alpha = 1, result_type beta = 1);
857
858 result_type a() const;
859 result_type b() const;
860
861 friend bool operator==(const param_type& x, const param_type& y);
862 friend bool operator!=(const param_type& x, const param_type& y);
863 };
864
865 // constructor and reset functions
866 explicit weibull_distribution(result_type a = 1, result_type b = 1);
867 explicit weibull_distribution(const param_type& parm);
868 void reset();
869
870 // generating functions
871 template<class URNG> result_type operator()(URNG& g);
872 template<class URNG> result_type operator()(URNG& g, const param_type& parm);
873
874 // property functions
875 result_type a() const;
876 result_type b() const;
877
878 param_type param() const;
879 void param(const param_type& parm);
880
881 result_type min() const;
882 result_type max() const;
883
884
885 friend bool operator==(const weibull_distribution& x,
886 const weibull_distribution& y);
887 friend bool operator!=(const weibull_distribution& x,
888 const weibull_distribution& y);
889
890 template <class charT, class traits>
891 friend
892 basic_ostream<charT, traits>&
893 operator<<(basic_ostream<charT, traits>& os,
894 const weibull_distribution& x);
895
896 template <class charT, class traits>
897 friend
898 basic_istream<charT, traits>&
899 operator>>(basic_istream<charT, traits>& is,
900 weibull_distribution& x);
901};
Howard Hinnantbc8d3f92010-05-11 19:42:16 +0000902
903template<class RealType = double>
904 class extreme_value_distribution;
905
906template<class RealType = double>
Howard Hinnanta64111c2010-05-12 21:02:31 +0000907class normal_distribution
908{
909public:
910 // types
911 typedef RealType result_type;
912
913 class param_type
914 {
915 public:
916 typedef normal_distribution distribution_type;
917
918 explicit param_type(result_type mean = 0, result_type stddev = 1);
919
920 result_type mean() const;
921 result_type stddev() const;
922
923 friend bool operator==(const param_type& x, const param_type& y);
924 friend bool operator!=(const param_type& x, const param_type& y);
925 };
926
927 // constructors and reset functions
928 explicit normal_distribution(result_type mean = 0, result_type stddev = 1);
929 explicit normal_distribution(const param_type& parm);
930 void reset();
931
932 // generating functions
933 template<class URNG> result_type operator()(URNG& g);
934 template<class URNG> result_type operator()(URNG& g, const param_type& parm);
935
936 // property functions
937 result_type mean() const;
938 result_type stddev() const;
939
940 param_type param() const;
941 void param(const param_type& parm);
942
943 result_type min() const;
944 result_type max() const;
945
946 friend bool operator==(const normal_distribution& x,
947 const normal_distribution& y);
948 friend bool operator!=(const normal_distribution& x,
949 const normal_distribution& y);
950
951 template <class charT, class traits>
952 friend
953 basic_ostream<charT, traits>&
954 operator<<(basic_ostream<charT, traits>& os,
955 const normal_distribution& x);
956
957 template <class charT, class traits>
958 friend
959 basic_istream<charT, traits>&
960 operator>>(basic_istream<charT, traits>& is,
961 normal_distribution& x);
962};
Howard Hinnantbc8d3f92010-05-11 19:42:16 +0000963
964template<class RealType = double>
965 class lognormal_distribution;
966
967template<class RealType = double>
Howard Hinnant97dc2f32010-05-15 23:36:00 +0000968class chi_squared_distribution
969{
970public:
971 // types
972 typedef RealType result_type;
973
974 class param_type
975 {
976 public:
977 typedef chi_squared_distribution distribution_type;
978
979 explicit param_type(result_type n = 1);
980
981 result_type n() const;
982
983 friend bool operator==(const param_type& x, const param_type& y);
984 friend bool operator!=(const param_type& x, const param_type& y);
985 };
986
987 // constructor and reset functions
988 explicit chi_squared_distribution(result_type n = 1);
989 explicit chi_squared_distribution(const param_type& parm);
990 void reset();
991
992 // generating functions
993 template<class URNG> result_type operator()(URNG& g);
994 template<class URNG> result_type operator()(URNG& g, const param_type& parm);
995
996 // property functions
997 result_type n() const;
998
999 param_type param() const;
1000 void param(const param_type& parm);
1001
1002 result_type min() const;
1003 result_type max() const;
1004
1005
1006 friend bool operator==(const chi_squared_distribution& x,
1007 const chi_squared_distribution& y);
1008 friend bool operator!=(const chi_squared_distribution& x,
1009 const chi_squared_distribution& y);
1010
1011 template <class charT, class traits>
1012 friend
1013 basic_ostream<charT, traits>&
1014 operator<<(basic_ostream<charT, traits>& os,
1015 const chi_squared_distribution& x);
1016
1017 template <class charT, class traits>
1018 friend
1019 basic_istream<charT, traits>&
1020 operator>>(basic_istream<charT, traits>& is,
1021 chi_squared_distribution& x);
1022};
Howard Hinnantbc8d3f92010-05-11 19:42:16 +00001023
1024template<class RealType = double>
1025 class cauchy_distribution;
1026
1027template<class RealType = double>
1028 class fisher_f_distribution;
1029
1030template<class RealType = double>
1031 class student_t_distribution;
1032
1033template<class IntType = int>
1034 class discrete_distribution;
1035
1036template<class RealType = double>
1037 class piecewise_constant_distribution;
1038
1039template<class RealType = double>
1040 class piecewise_linear_distribution;
1041
1042} // std
1043*/
1044
1045#include <__config>
1046#include <cstddef>
1047#include <type_traits>
1048#include <initializer_list>
1049#include <cstdint>
1050#include <limits>
1051#include <algorithm>
1052#include <vector>
1053#include <string>
1054#include <istream>
1055#include <ostream>
Howard Hinnant30a840f2010-05-12 17:08:57 +00001056#include <cmath>
Howard Hinnantbc8d3f92010-05-11 19:42:16 +00001057
1058#pragma GCC system_header
1059
1060_LIBCPP_BEGIN_NAMESPACE_STD
1061
1062// linear_congruential_engine
1063
1064template <unsigned long long __a, unsigned long long __c,
1065 unsigned long long __m, unsigned long long _M,
1066 bool _MightOverflow = (__a != 0 && __m != 0 && __m-1 > (_M-__c)/__a)>
1067struct __lce_ta;
1068
1069// 64
1070
1071template <unsigned long long __a, unsigned long long __c, unsigned long long __m>
1072struct __lce_ta<__a, __c, __m, (unsigned long long)(~0), true>
1073{
1074 typedef unsigned long long result_type;
1075 static result_type next(result_type __x)
1076 {
1077 // Schrage's algorithm
1078 const result_type __q = __m / __a;
1079 const result_type __r = __m % __a;
1080 const result_type __t0 = __a * (__x % __q);
1081 const result_type __t1 = __r * (__x / __q);
1082 __x = __t0 + (__t0 < __t1) * __m - __t1;
1083 __x += __c - (__x >= __m - __c) * __m;
1084 return __x;
1085 }
1086};
1087
1088template <unsigned long long __a, unsigned long long __m>
1089struct __lce_ta<__a, 0, __m, (unsigned long long)(~0), true>
1090{
1091 typedef unsigned long long result_type;
1092 static result_type next(result_type __x)
1093 {
1094 // Schrage's algorithm
1095 const result_type __q = __m / __a;
1096 const result_type __r = __m % __a;
1097 const result_type __t0 = __a * (__x % __q);
1098 const result_type __t1 = __r * (__x / __q);
1099 __x = __t0 + (__t0 < __t1) * __m - __t1;
1100 return __x;
1101 }
1102};
1103
1104template <unsigned long long __a, unsigned long long __c, unsigned long long __m>
1105struct __lce_ta<__a, __c, __m, (unsigned long long)(~0), false>
1106{
1107 typedef unsigned long long result_type;
1108 static result_type next(result_type __x)
1109 {
1110 return (__a * __x + __c) % __m;
1111 }
1112};
1113
1114template <unsigned long long __a, unsigned long long __c>
1115struct __lce_ta<__a, __c, 0, (unsigned long long)(~0), false>
1116{
1117 typedef unsigned long long result_type;
1118 static result_type next(result_type __x)
1119 {
1120 return __a * __x + __c;
1121 }
1122};
1123
1124// 32
1125
1126template <unsigned long long _A, unsigned long long _C, unsigned long long _M>
1127struct __lce_ta<_A, _C, _M, unsigned(~0), true>
1128{
1129 typedef unsigned result_type;
1130 static result_type next(result_type __x)
1131 {
1132 const result_type __a = static_cast<result_type>(_A);
1133 const result_type __c = static_cast<result_type>(_C);
1134 const result_type __m = static_cast<result_type>(_M);
1135 // Schrage's algorithm
1136 const result_type __q = __m / __a;
1137 const result_type __r = __m % __a;
1138 const result_type __t0 = __a * (__x % __q);
1139 const result_type __t1 = __r * (__x / __q);
1140 __x = __t0 + (__t0 < __t1) * __m - __t1;
1141 __x += __c - (__x >= __m - __c) * __m;
1142 return __x;
1143 }
1144};
1145
1146template <unsigned long long _A, unsigned long long _M>
1147struct __lce_ta<_A, 0, _M, unsigned(~0), true>
1148{
1149 typedef unsigned result_type;
1150 static result_type next(result_type __x)
1151 {
1152 const result_type __a = static_cast<result_type>(_A);
1153 const result_type __m = static_cast<result_type>(_M);
1154 // Schrage's algorithm
1155 const result_type __q = __m / __a;
1156 const result_type __r = __m % __a;
1157 const result_type __t0 = __a * (__x % __q);
1158 const result_type __t1 = __r * (__x / __q);
1159 __x = __t0 + (__t0 < __t1) * __m - __t1;
1160 return __x;
1161 }
1162};
1163
1164template <unsigned long long _A, unsigned long long _C, unsigned long long _M>
1165struct __lce_ta<_A, _C, _M, unsigned(~0), false>
1166{
1167 typedef unsigned result_type;
1168 static result_type next(result_type __x)
1169 {
1170 const result_type __a = static_cast<result_type>(_A);
1171 const result_type __c = static_cast<result_type>(_C);
1172 const result_type __m = static_cast<result_type>(_M);
1173 return (__a * __x + __c) % __m;
1174 }
1175};
1176
1177template <unsigned long long _A, unsigned long long _C>
1178struct __lce_ta<_A, _C, 0, unsigned(~0), false>
1179{
1180 typedef unsigned result_type;
1181 static result_type next(result_type __x)
1182 {
1183 const result_type __a = static_cast<result_type>(_A);
1184 const result_type __c = static_cast<result_type>(_C);
1185 return __a * __x + __c;
1186 }
1187};
1188
1189// 16
1190
1191template <unsigned long long __a, unsigned long long __c, unsigned long long __m, bool __b>
1192struct __lce_ta<__a, __c, __m, (unsigned short)(~0), __b>
1193{
1194 typedef unsigned short result_type;
1195 static result_type next(result_type __x)
1196 {
1197 return static_cast<result_type>(__lce_ta<__a, __c, __m, unsigned(~0)>::next(__x));
1198 }
1199};
1200
1201template <class _UIntType, _UIntType __a, _UIntType __c, _UIntType __m>
1202class linear_congruential_engine;
1203
1204template <class _CharT, class _Traits,
1205 class _U, _U _A, _U _C, _U _N>
1206basic_ostream<_CharT, _Traits>&
1207operator<<(basic_ostream<_CharT, _Traits>& __os,
1208 const linear_congruential_engine<_U, _A, _C, _N>&);
1209
1210template <class _CharT, class _Traits,
1211 class _U, _U _A, _U _C, _U _N>
1212basic_istream<_CharT, _Traits>&
1213operator>>(basic_istream<_CharT, _Traits>& __is,
1214 linear_congruential_engine<_U, _A, _C, _N>& __x);
1215
1216template <class _UIntType, _UIntType __a, _UIntType __c, _UIntType __m>
1217class linear_congruential_engine
1218{
1219public:
1220 // types
1221 typedef _UIntType result_type;
1222
1223private:
1224 result_type __x_;
1225
1226 static const result_type _M = result_type(~0);
1227
1228 static_assert(__m == 0 || __a < __m, "linear_congruential_engine invalid parameters");
1229 static_assert(__m == 0 || __c < __m, "linear_congruential_engine invalid parameters");
1230public:
1231 static const result_type _Min = __c == 0u ? 1u: 0u;
1232 static const result_type _Max = __m - 1u;
1233 static_assert(_Min < _Max, "linear_congruential_engine invalid parameters");
1234
1235 // engine characteristics
1236 static const/*expr*/ result_type multiplier = __a;
1237 static const/*expr*/ result_type increment = __c;
1238 static const/*expr*/ result_type modulus = __m;
1239 static const/*expr*/ result_type min() {return _Min;}
1240 static const/*expr*/ result_type max() {return _Max;}
1241 static const/*expr*/ result_type default_seed = 1u;
1242
1243 // constructors and seeding functions
1244 explicit linear_congruential_engine(result_type __s = default_seed)
1245 {seed(__s);}
1246 template<class _Sseq> explicit linear_congruential_engine(_Sseq& __q)
1247 {seed(__q);}
1248 void seed(result_type __s = default_seed)
1249 {seed(integral_constant<bool, __m == 0>(),
1250 integral_constant<bool, __c == 0>(), __s);}
1251 template<class _Sseq>
1252 typename enable_if
1253 <
1254 !is_convertible<_Sseq, result_type>::value,
1255 void
1256 >::type
1257 seed(_Sseq& __q)
1258 {__seed(__q, integral_constant<unsigned,
1259 1 + (__m == 0 ? (sizeof(result_type) * __CHAR_BIT__ - 1)/32
1260 : (__m-1) / 0x100000000ull)>());}
1261
1262 // generating functions
1263 result_type operator()()
1264 {return __x_ = static_cast<result_type>(__lce_ta<__a, __c, __m, _M>::next(__x_));}
1265 void discard(unsigned long long __z) {for (; __z; --__z) operator()();}
1266
1267 friend bool operator==(const linear_congruential_engine& __x,
1268 const linear_congruential_engine& __y)
1269 {return __x.__x_ == __y.__x_;}
1270 friend bool operator!=(const linear_congruential_engine& __x,
1271 const linear_congruential_engine& __y)
1272 {return !(__x == __y);}
1273
1274private:
1275
1276 void seed(true_type, true_type, result_type __s) {__x_ = __s == 0 ? 1 : __s;}
1277 void seed(true_type, false_type, result_type __s) {__x_ = __s;}
1278 void seed(false_type, true_type, result_type __s) {__x_ = __s % __m == 0 ?
1279 1 : __s % __m;}
1280 void seed(false_type, false_type, result_type __s) {__x_ = __s % __m;}
1281
1282 template<class _Sseq>
1283 void __seed(_Sseq& __q, integral_constant<unsigned, 1>);
1284 template<class _Sseq>
1285 void __seed(_Sseq& __q, integral_constant<unsigned, 2>);
1286
1287 template <class _CharT, class _Traits,
1288 class _U, _U _A, _U _C, _U _N>
1289 friend
1290 basic_ostream<_CharT, _Traits>&
1291 operator<<(basic_ostream<_CharT, _Traits>& __os,
1292 const linear_congruential_engine<_U, _A, _C, _N>&);
1293
1294 template <class _CharT, class _Traits,
1295 class _U, _U _A, _U _C, _U _N>
1296 friend
1297 basic_istream<_CharT, _Traits>&
1298 operator>>(basic_istream<_CharT, _Traits>& __is,
1299 linear_congruential_engine<_U, _A, _C, _N>& __x);
1300};
1301
1302template <class _UIntType, _UIntType __a, _UIntType __c, _UIntType __m>
1303template<class _Sseq>
1304void
1305linear_congruential_engine<_UIntType, __a, __c, __m>::__seed(_Sseq& __q,
1306 integral_constant<unsigned, 1>)
1307{
1308 const unsigned __k = 1;
1309 uint32_t __ar[__k+3];
1310 __q.generate(__ar, __ar + __k + 3);
1311 result_type __s = static_cast<result_type>(__ar[3] % __m);
1312 __x_ = __c == 0 && __s == 0 ? result_type(1) : __s;
1313}
1314
1315template <class _UIntType, _UIntType __a, _UIntType __c, _UIntType __m>
1316template<class _Sseq>
1317void
1318linear_congruential_engine<_UIntType, __a, __c, __m>::__seed(_Sseq& __q,
1319 integral_constant<unsigned, 2>)
1320{
1321 const unsigned __k = 2;
1322 uint32_t __ar[__k+3];
1323 __q.generate(__ar, __ar + __k + 3);
1324 result_type __s = static_cast<result_type>((__ar[3] +
1325 (uint64_t)__ar[4] << 32) % __m);
1326 __x_ = __c == 0 && __s == 0 ? result_type(1) : __s;
1327}
1328
1329template <class _CharT, class _Traits>
1330class __save_flags
1331{
1332 typedef basic_ios<_CharT, _Traits> __stream_type;
1333 typedef typename __stream_type::fmtflags fmtflags;
1334
1335 __stream_type& __stream_;
1336 fmtflags __fmtflags_;
1337 _CharT __fill_;
1338
1339 __save_flags(const __save_flags&);
1340 __save_flags& operator=(const __save_flags&);
1341public:
1342 explicit __save_flags(__stream_type& __stream)
1343 : __stream_(__stream),
1344 __fmtflags_(__stream.flags()),
1345 __fill_(__stream.fill())
1346 {}
1347 ~__save_flags()
1348 {
1349 __stream_.flags(__fmtflags_);
1350 __stream_.fill(__fill_);
1351 }
1352};
1353
1354template <class _CharT, class _Traits,
1355 class _UIntType, _UIntType __a, _UIntType __c, _UIntType __m>
1356inline
1357basic_ostream<_CharT, _Traits>&
1358operator<<(basic_ostream<_CharT, _Traits>& __os,
1359 const linear_congruential_engine<_UIntType, __a, __c, __m>& __x)
1360{
1361 __save_flags<_CharT, _Traits> _(__os);
1362 __os.flags(ios_base::dec | ios_base::left);
1363 __os.fill(__os.widen(' '));
1364 return __os << __x.__x_;
1365}
1366
1367template <class _CharT, class _Traits,
1368 class _UIntType, _UIntType __a, _UIntType __c, _UIntType __m>
1369basic_istream<_CharT, _Traits>&
1370operator>>(basic_istream<_CharT, _Traits>& __is,
1371 linear_congruential_engine<_UIntType, __a, __c, __m>& __x)
1372{
1373 __save_flags<_CharT, _Traits> _(__is);
1374 __is.flags(ios_base::dec | ios_base::skipws);
1375 _UIntType __t;
1376 __is >> __t;
1377 if (!__is.fail())
1378 __x.__x_ = __t;
1379 return __is;
1380}
1381
1382typedef linear_congruential_engine<uint_fast32_t, 16807, 0, 2147483647>
1383 minstd_rand0;
1384typedef minstd_rand0 default_random_engine;
1385typedef linear_congruential_engine<uint_fast32_t, 48271, 0, 2147483647>
1386 minstd_rand;
1387// mersenne_twister_engine
1388
1389template <class _UIntType, size_t __w, size_t __n, size_t __m, size_t __r,
1390 _UIntType __a, size_t __u, _UIntType __d, size_t __s,
1391 _UIntType __b, size_t __t, _UIntType __c, size_t __l, _UIntType __f>
1392class mersenne_twister_engine;
1393
1394template <class _UI, size_t _W, size_t _N, size_t _M, size_t _R,
1395 _UI _A, size_t _U, _UI _D, size_t _S,
1396 _UI _B, size_t _T, _UI _C, size_t _L, _UI _F>
1397bool
1398operator==(const mersenne_twister_engine<_UI, _W, _N, _M, _R, _A, _U, _D, _S,
1399 _B, _T, _C, _L, _F>& __x,
1400 const mersenne_twister_engine<_UI, _W, _N, _M, _R, _A, _U, _D, _S,
1401 _B, _T, _C, _L, _F>& __y);
1402
1403template <class _UI, size_t _W, size_t _N, size_t _M, size_t _R,
1404 _UI _A, size_t _U, _UI _D, size_t _S,
1405 _UI _B, size_t _T, _UI _C, size_t _L, _UI _F>
1406bool
1407operator!=(const mersenne_twister_engine<_UI, _W, _N, _M, _R, _A, _U, _D, _S,
1408 _B, _T, _C, _L, _F>& __x,
1409 const mersenne_twister_engine<_UI, _W, _N, _M, _R, _A, _U, _D, _S,
1410 _B, _T, _C, _L, _F>& __y);
1411
1412template <class _CharT, class _Traits,
1413 class _UI, size_t _W, size_t _N, size_t _M, size_t _R,
1414 _UI _A, size_t _U, _UI _D, size_t _S,
1415 _UI _B, size_t _T, _UI _C, size_t _L, _UI _F>
1416basic_ostream<_CharT, _Traits>&
1417operator<<(basic_ostream<_CharT, _Traits>& __os,
1418 const mersenne_twister_engine<_UI, _W, _N, _M, _R, _A, _U, _D, _S,
1419 _B, _T, _C, _L, _F>& __x);
1420
1421template <class _CharT, class _Traits,
1422 class _UI, size_t _W, size_t _N, size_t _M, size_t _R,
1423 _UI _A, size_t _U, _UI _D, size_t _S,
1424 _UI _B, size_t _T, _UI _C, size_t _L, _UI _F>
1425basic_istream<_CharT, _Traits>&
1426operator>>(basic_istream<_CharT, _Traits>& __is,
1427 mersenne_twister_engine<_UI, _W, _N, _M, _R, _A, _U, _D, _S,
1428 _B, _T, _C, _L, _F>& __x);
1429
1430template <class _UIntType, size_t __w, size_t __n, size_t __m, size_t __r,
1431 _UIntType __a, size_t __u, _UIntType __d, size_t __s,
1432 _UIntType __b, size_t __t, _UIntType __c, size_t __l, _UIntType __f>
1433class mersenne_twister_engine
1434{
1435public:
1436 // types
1437 typedef _UIntType result_type;
1438
1439private:
1440 result_type __x_[__n];
1441 size_t __i_;
1442
1443 static_assert( 0 < __m, "mersenne_twister_engine invalid parameters");
1444 static_assert(__m <= __n, "mersenne_twister_engine invalid parameters");
1445 static const result_type _Dt = numeric_limits<result_type>::digits;
1446 static_assert(__w <= _Dt, "mersenne_twister_engine invalid parameters");
1447 static_assert( 2 <= __w, "mersenne_twister_engine invalid parameters");
1448 static_assert(__r <= __w, "mersenne_twister_engine invalid parameters");
1449 static_assert(__u <= __w, "mersenne_twister_engine invalid parameters");
1450 static_assert(__s <= __w, "mersenne_twister_engine invalid parameters");
1451 static_assert(__t <= __w, "mersenne_twister_engine invalid parameters");
1452 static_assert(__l <= __w, "mersenne_twister_engine invalid parameters");
1453public:
1454 static const result_type _Min = 0;
1455 static const result_type _Max = __w == _Dt ? result_type(~0) :
1456 (result_type(1) << __w) - result_type(1);
1457 static_assert(_Min < _Max, "mersenne_twister_engine invalid parameters");
1458 static_assert(__a <= _Max, "mersenne_twister_engine invalid parameters");
1459 static_assert(__b <= _Max, "mersenne_twister_engine invalid parameters");
1460 static_assert(__c <= _Max, "mersenne_twister_engine invalid parameters");
1461 static_assert(__d <= _Max, "mersenne_twister_engine invalid parameters");
1462 static_assert(__f <= _Max, "mersenne_twister_engine invalid parameters");
1463
1464 // engine characteristics
1465 static const/*expr*/ size_t word_size = __w;
1466 static const/*expr*/ size_t state_size = __n;
1467 static const/*expr*/ size_t shift_size = __m;
1468 static const/*expr*/ size_t mask_bits = __r;
1469 static const/*expr*/ result_type xor_mask = __a;
1470 static const/*expr*/ size_t tempering_u = __u;
1471 static const/*expr*/ result_type tempering_d = __d;
1472 static const/*expr*/ size_t tempering_s = __s;
1473 static const/*expr*/ result_type tempering_b = __b;
1474 static const/*expr*/ size_t tempering_t = __t;
1475 static const/*expr*/ result_type tempering_c = __c;
1476 static const/*expr*/ size_t tempering_l = __l;
1477 static const/*expr*/ result_type initialization_multiplier = __f;
1478 static const/*expr*/ result_type min() { return _Min; }
1479 static const/*expr*/ result_type max() { return _Max; }
1480 static const/*expr*/ result_type default_seed = 5489u;
1481
1482 // constructors and seeding functions
1483 explicit mersenne_twister_engine(result_type __sd = default_seed)
1484 {seed(__sd);}
1485 template<class _Sseq> explicit mersenne_twister_engine(_Sseq& __q)
1486 {seed(__q);}
1487 void seed(result_type __sd = default_seed);
1488 template<class _Sseq>
1489 typename enable_if
1490 <
1491 !is_convertible<_Sseq, result_type>::value,
1492 void
1493 >::type
1494 seed(_Sseq& __q)
1495 {__seed(__q, integral_constant<unsigned, 1 + (__w - 1) / 32>());}
1496
1497 // generating functions
1498 result_type operator()();
1499 void discard(unsigned long long __z) {for (; __z; --__z) operator()();}
1500
1501 template <class _UI, size_t _W, size_t _N, size_t _M, size_t _R,
1502 _UI _A, size_t _U, _UI _D, size_t _S,
1503 _UI _B, size_t _T, _UI _C, size_t _L, _UI _F>
1504 friend
1505 bool
1506 operator==(const mersenne_twister_engine<_UI, _W, _N, _M, _R, _A, _U, _D, _S,
1507 _B, _T, _C, _L, _F>& __x,
1508 const mersenne_twister_engine<_UI, _W, _N, _M, _R, _A, _U, _D, _S,
1509 _B, _T, _C, _L, _F>& __y);
1510
1511 template <class _UI, size_t _W, size_t _N, size_t _M, size_t _R,
1512 _UI _A, size_t _U, _UI _D, size_t _S,
1513 _UI _B, size_t _T, _UI _C, size_t _L, _UI _F>
1514 friend
1515 bool
1516 operator!=(const mersenne_twister_engine<_UI, _W, _N, _M, _R, _A, _U, _D, _S,
1517 _B, _T, _C, _L, _F>& __x,
1518 const mersenne_twister_engine<_UI, _W, _N, _M, _R, _A, _U, _D, _S,
1519 _B, _T, _C, _L, _F>& __y);
1520
1521 template <class _CharT, class _Traits,
1522 class _UI, size_t _W, size_t _N, size_t _M, size_t _R,
1523 _UI _A, size_t _U, _UI _D, size_t _S,
1524 _UI _B, size_t _T, _UI _C, size_t _L, _UI _F>
1525 friend
1526 basic_ostream<_CharT, _Traits>&
1527 operator<<(basic_ostream<_CharT, _Traits>& __os,
1528 const mersenne_twister_engine<_UI, _W, _N, _M, _R, _A, _U, _D, _S,
1529 _B, _T, _C, _L, _F>& __x);
1530
1531 template <class _CharT, class _Traits,
1532 class _UI, size_t _W, size_t _N, size_t _M, size_t _R,
1533 _UI _A, size_t _U, _UI _D, size_t _S,
1534 _UI _B, size_t _T, _UI _C, size_t _L, _UI _F>
1535 friend
1536 basic_istream<_CharT, _Traits>&
1537 operator>>(basic_istream<_CharT, _Traits>& __is,
1538 mersenne_twister_engine<_UI, _W, _N, _M, _R, _A, _U, _D, _S,
1539 _B, _T, _C, _L, _F>& __x);
1540private:
1541
1542 template<class _Sseq>
1543 void __seed(_Sseq& __q, integral_constant<unsigned, 1>);
1544 template<class _Sseq>
1545 void __seed(_Sseq& __q, integral_constant<unsigned, 2>);
1546
1547 template <size_t __count>
1548 static
1549 typename enable_if
1550 <
1551 __count < __w,
1552 result_type
1553 >::type
1554 __lshift(result_type __x) {return (__x << __count) & _Max;}
1555
1556 template <size_t __count>
1557 static
1558 typename enable_if
1559 <
1560 (__count >= __w),
1561 result_type
1562 >::type
1563 __lshift(result_type __x) {return result_type(0);}
1564
1565 template <size_t __count>
1566 static
1567 typename enable_if
1568 <
1569 __count < _Dt,
1570 result_type
1571 >::type
1572 __rshift(result_type __x) {return __x >> __count;}
1573
1574 template <size_t __count>
1575 static
1576 typename enable_if
1577 <
1578 (__count >= _Dt),
1579 result_type
1580 >::type
1581 __rshift(result_type __x) {return result_type(0);}
1582};
1583
1584template <class _UIntType, size_t __w, size_t __n, size_t __m, size_t __r,
1585 _UIntType __a, size_t __u, _UIntType __d, size_t __s,
1586 _UIntType __b, size_t __t, _UIntType __c, size_t __l, _UIntType __f>
1587void
1588mersenne_twister_engine<_UIntType, __w, __n, __m, __r, __a, __u, __d, __s, __b,
1589 __t, __c, __l, __f>::seed(result_type __sd)
1590{ // __w >= 2
1591 __x_[0] = __sd & _Max;
1592 for (size_t __i = 1; __i < __n; ++__i)
1593 __x_[__i] = (__f * (__x_[__i-1] ^ __rshift<__w - 2>(__x_[__i-1])) + __i) & _Max;
1594 __i_ = 0;
1595}
1596
1597template <class _UIntType, size_t __w, size_t __n, size_t __m, size_t __r,
1598 _UIntType __a, size_t __u, _UIntType __d, size_t __s,
1599 _UIntType __b, size_t __t, _UIntType __c, size_t __l, _UIntType __f>
1600template<class _Sseq>
1601void
1602mersenne_twister_engine<_UIntType, __w, __n, __m, __r, __a, __u, __d, __s, __b,
1603 __t, __c, __l, __f>::__seed(_Sseq& __q, integral_constant<unsigned, 1>)
1604{
1605 const unsigned __k = 1;
1606 uint32_t __ar[__n * __k];
1607 __q.generate(__ar, __ar + __n * __k);
1608 for (size_t __i = 0; __i < __n; ++__i)
1609 __x_[__i] = static_cast<result_type>(__ar[__i] & _Max);
1610 const result_type __mask = __r == _Dt ? result_type(~0) :
1611 (result_type(1) << __r) - result_type(1);
1612 __i_ = 0;
1613 if ((__x_[0] & ~__mask) == 0)
1614 {
1615 for (size_t __i = 1; __i < __n; ++__i)
1616 if (__x_[__i] != 0)
1617 return;
1618 __x_[0] = _Max;
1619 }
1620}
1621
1622template <class _UIntType, size_t __w, size_t __n, size_t __m, size_t __r,
1623 _UIntType __a, size_t __u, _UIntType __d, size_t __s,
1624 _UIntType __b, size_t __t, _UIntType __c, size_t __l, _UIntType __f>
1625template<class _Sseq>
1626void
1627mersenne_twister_engine<_UIntType, __w, __n, __m, __r, __a, __u, __d, __s, __b,
1628 __t, __c, __l, __f>::__seed(_Sseq& __q, integral_constant<unsigned, 2>)
1629{
1630 const unsigned __k = 2;
1631 uint32_t __ar[__n * __k];
1632 __q.generate(__ar, __ar + __n * __k);
1633 for (size_t __i = 0; __i < __n; ++__i)
1634 __x_[__i] = static_cast<result_type>(
1635 (__ar[2 * __i] + ((uint64_t)__ar[2 * __i + 1] << 32)) & _Max);
1636 const result_type __mask = __r == _Dt ? result_type(~0) :
1637 (result_type(1) << __r) - result_type(1);
1638 __i_ = 0;
1639 if ((__x_[0] & ~__mask) == 0)
1640 {
1641 for (size_t __i = 1; __i < __n; ++__i)
1642 if (__x_[__i] != 0)
1643 return;
1644 __x_[0] = _Max;
1645 }
1646}
1647
1648template <class _UIntType, size_t __w, size_t __n, size_t __m, size_t __r,
1649 _UIntType __a, size_t __u, _UIntType __d, size_t __s,
1650 _UIntType __b, size_t __t, _UIntType __c, size_t __l, _UIntType __f>
1651_UIntType
1652mersenne_twister_engine<_UIntType, __w, __n, __m, __r, __a, __u, __d, __s, __b,
1653 __t, __c, __l, __f>::operator()()
1654{
1655 const size_t __j = (__i_ + 1) % __n;
1656 const result_type __mask = __r == _Dt ? result_type(~0) :
1657 (result_type(1) << __r) - result_type(1);
1658 const result_type _Y = (__x_[__i_] & ~__mask) | (__x_[__j] & __mask);
1659 const size_t __k = (__i_ + __m) % __n;
1660 __x_[__i_] = __x_[__k] ^ __rshift<1>(_Y) ^ (__a * (_Y & 1));
1661 result_type __z = __x_[__i_] ^ (__rshift<__u>(__x_[__i_]) & __d);
1662 __i_ = __j;
1663 __z ^= __lshift<__s>(__z) & __b;
1664 __z ^= __lshift<__t>(__z) & __c;
1665 return __z ^ __rshift<__l>(__z);
1666}
1667
1668template <class _UI, size_t _W, size_t _N, size_t _M, size_t _R,
1669 _UI _A, size_t _U, _UI _D, size_t _S,
1670 _UI _B, size_t _T, _UI _C, size_t _L, _UI _F>
1671bool
1672operator==(const mersenne_twister_engine<_UI, _W, _N, _M, _R, _A, _U, _D, _S,
1673 _B, _T, _C, _L, _F>& __x,
1674 const mersenne_twister_engine<_UI, _W, _N, _M, _R, _A, _U, _D, _S,
1675 _B, _T, _C, _L, _F>& __y)
1676{
1677 if (__x.__i_ == __y.__i_)
1678 return _STD::equal(__x.__x_, __x.__x_ + _N, __y.__x_);
1679 if (__x.__i_ == 0 || __y.__i_ == 0)
1680 {
1681 size_t __j = _STD::min(_N - __x.__i_, _N - __y.__i_);
1682 if (!_STD::equal(__x.__x_ + __x.__i_, __x.__x_ + __x.__i_ + __j,
1683 __y.__x_ + __y.__i_))
1684 return false;
1685 if (__x.__i_ == 0)
1686 return _STD::equal(__x.__x_ + __j, __x.__x_ + _N, __y.__x_);
1687 return _STD::equal(__x.__x_, __x.__x_ + (_N - __j), __y.__x_ + __j);
1688 }
1689 if (__x.__i_ < __y.__i_)
1690 {
1691 size_t __j = _N - __y.__i_;
1692 if (!_STD::equal(__x.__x_ + __x.__i_, __x.__x_ + (__x.__i_ + __j),
1693 __y.__x_ + __y.__i_))
1694 return false;
1695 if (!_STD::equal(__x.__x_ + (__x.__i_ + __j), __x.__x_ + _N,
1696 __y.__x_))
1697 return false;
1698 return _STD::equal(__x.__x_, __x.__x_ + __x.__i_,
1699 __y.__x_ + (_N - (__x.__i_ + __j)));
1700 }
1701 size_t __j = _N - __x.__i_;
1702 if (!_STD::equal(__y.__x_ + __y.__i_, __y.__x_ + (__y.__i_ + __j),
1703 __x.__x_ + __x.__i_))
1704 return false;
1705 if (!_STD::equal(__y.__x_ + (__y.__i_ + __j), __y.__x_ + _N,
1706 __x.__x_))
1707 return false;
1708 return _STD::equal(__y.__x_, __y.__x_ + __y.__i_,
1709 __x.__x_ + (_N - (__y.__i_ + __j)));
1710}
1711
1712template <class _UI, size_t _W, size_t _N, size_t _M, size_t _R,
1713 _UI _A, size_t _U, _UI _D, size_t _S,
1714 _UI _B, size_t _T, _UI _C, size_t _L, _UI _F>
1715inline
1716bool
1717operator!=(const mersenne_twister_engine<_UI, _W, _N, _M, _R, _A, _U, _D, _S,
1718 _B, _T, _C, _L, _F>& __x,
1719 const mersenne_twister_engine<_UI, _W, _N, _M, _R, _A, _U, _D, _S,
1720 _B, _T, _C, _L, _F>& __y)
1721{
1722 return !(__x == __y);
1723}
1724
1725template <class _CharT, class _Traits,
1726 class _UI, size_t _W, size_t _N, size_t _M, size_t _R,
1727 _UI _A, size_t _U, _UI _D, size_t _S,
1728 _UI _B, size_t _T, _UI _C, size_t _L, _UI _F>
1729basic_ostream<_CharT, _Traits>&
1730operator<<(basic_ostream<_CharT, _Traits>& __os,
1731 const mersenne_twister_engine<_UI, _W, _N, _M, _R, _A, _U, _D, _S,
1732 _B, _T, _C, _L, _F>& __x)
1733{
1734 __save_flags<_CharT, _Traits> _(__os);
1735 __os.flags(ios_base::dec | ios_base::left);
1736 _CharT __sp = __os.widen(' ');
1737 __os.fill(__sp);
1738 __os << __x.__x_[__x.__i_];
1739 for (size_t __j = __x.__i_ + 1; __j < _N; ++__j)
1740 __os << __sp << __x.__x_[__j];
1741 for (size_t __j = 0; __j < __x.__i_; ++__j)
1742 __os << __sp << __x.__x_[__j];
1743 return __os;
1744}
1745
1746template <class _CharT, class _Traits,
1747 class _UI, size_t _W, size_t _N, size_t _M, size_t _R,
1748 _UI _A, size_t _U, _UI _D, size_t _S,
1749 _UI _B, size_t _T, _UI _C, size_t _L, _UI _F>
1750basic_istream<_CharT, _Traits>&
1751operator>>(basic_istream<_CharT, _Traits>& __is,
1752 mersenne_twister_engine<_UI, _W, _N, _M, _R, _A, _U, _D, _S,
1753 _B, _T, _C, _L, _F>& __x)
1754{
1755 __save_flags<_CharT, _Traits> _(__is);
1756 __is.flags(ios_base::dec | ios_base::skipws);
1757 _UI __t[_N];
1758 for (size_t __i = 0; __i < _N; ++__i)
1759 __is >> __t[__i];
1760 if (!__is.fail())
1761 {
1762 for (size_t __i = 0; __i < _N; ++__i)
1763 __x.__x_[__i] = __t[__i];
1764 __x.__i_ = 0;
1765 }
1766 return __is;
1767}
1768
1769typedef mersenne_twister_engine<uint_fast32_t, 32, 624, 397, 31,
1770 0x9908b0df, 11, 0xffffffff,
1771 7, 0x9d2c5680,
1772 15, 0xefc60000,
1773 18, 1812433253> mt19937;
1774typedef mersenne_twister_engine<uint_fast64_t, 64, 312, 156, 31,
1775 0xb5026f5aa96619e9ULL, 29, 0x5555555555555555ULL,
1776 17, 0x71d67fffeda60000ULL,
1777 37, 0xfff7eee000000000ULL,
1778 43, 6364136223846793005ULL> mt19937_64;
1779
1780// subtract_with_carry_engine
1781
1782template<class _UIntType, size_t __w, size_t __s, size_t __r>
1783class subtract_with_carry_engine;
1784
1785template<class _UI, size_t _W, size_t _S, size_t _R>
1786bool
1787operator==(
1788 const subtract_with_carry_engine<_UI, _W, _S, _R>& __x,
1789 const subtract_with_carry_engine<_UI, _W, _S, _R>& __y);
1790
1791template<class _UI, size_t _W, size_t _S, size_t _R>
1792bool
1793operator!=(
1794 const subtract_with_carry_engine<_UI, _W, _S, _R>& __x,
1795 const subtract_with_carry_engine<_UI, _W, _S, _R>& __y);
1796
1797template <class _CharT, class _Traits,
1798 class _UI, size_t _W, size_t _S, size_t _R>
1799basic_ostream<_CharT, _Traits>&
1800operator<<(basic_ostream<_CharT, _Traits>& __os,
1801 const subtract_with_carry_engine<_UI, _W, _S, _R>& __x);
1802
1803template <class _CharT, class _Traits,
1804 class _UI, size_t _W, size_t _S, size_t _R>
1805basic_istream<_CharT, _Traits>&
1806operator>>(basic_istream<_CharT, _Traits>& __is,
1807 subtract_with_carry_engine<_UI, _W, _S, _R>& __x);
1808
1809template<class _UIntType, size_t __w, size_t __s, size_t __r>
1810class subtract_with_carry_engine
1811{
1812public:
1813 // types
1814 typedef _UIntType result_type;
1815
1816private:
1817 result_type __x_[__r];
1818 result_type __c_;
1819 size_t __i_;
1820
1821 static const result_type _Dt = numeric_limits<result_type>::digits;
1822 static_assert( 0 < __w, "subtract_with_carry_engine invalid parameters");
1823 static_assert(__w <= _Dt, "subtract_with_carry_engine invalid parameters");
1824 static_assert( 0 < __s, "subtract_with_carry_engine invalid parameters");
1825 static_assert(__s < __r, "subtract_with_carry_engine invalid parameters");
1826public:
1827 static const result_type _Min = 0;
1828 static const result_type _Max = __w == _Dt ? result_type(~0) :
1829 (result_type(1) << __w) - result_type(1);
1830 static_assert(_Min < _Max, "subtract_with_carry_engine invalid parameters");
1831
1832 // engine characteristics
1833 static const/*expr*/ size_t word_size = __w;
1834 static const/*expr*/ size_t short_lag = __s;
1835 static const/*expr*/ size_t long_lag = __r;
1836 static const/*expr*/ result_type min() { return _Min; }
1837 static const/*expr*/ result_type max() { return _Max; }
1838 static const/*expr*/ result_type default_seed = 19780503u;
1839
1840 // constructors and seeding functions
1841 explicit subtract_with_carry_engine(result_type __sd = default_seed)
1842 {seed(__sd);}
1843 template<class _Sseq> explicit subtract_with_carry_engine(_Sseq& __q)
1844 {seed(__q);}
1845 void seed(result_type __sd = default_seed)
1846 {seed(__sd, integral_constant<unsigned, 1 + (__w - 1) / 32>());}
1847 template<class _Sseq>
1848 typename enable_if
1849 <
1850 !is_convertible<_Sseq, result_type>::value,
1851 void
1852 >::type
1853 seed(_Sseq& __q)
1854 {__seed(__q, integral_constant<unsigned, 1 + (__w - 1) / 32>());}
1855
1856 // generating functions
1857 result_type operator()();
1858 void discard(unsigned long long __z) {for (; __z; --__z) operator()();}
1859
1860 template<class _UI, size_t _W, size_t _S, size_t _R>
1861 friend
1862 bool
1863 operator==(
1864 const subtract_with_carry_engine<_UI, _W, _S, _R>& __x,
1865 const subtract_with_carry_engine<_UI, _W, _S, _R>& __y);
1866
1867 template<class _UI, size_t _W, size_t _S, size_t _R>
1868 friend
1869 bool
1870 operator!=(
1871 const subtract_with_carry_engine<_UI, _W, _S, _R>& __x,
1872 const subtract_with_carry_engine<_UI, _W, _S, _R>& __y);
1873
1874 template <class _CharT, class _Traits,
1875 class _UI, size_t _W, size_t _S, size_t _R>
1876 friend
1877 basic_ostream<_CharT, _Traits>&
1878 operator<<(basic_ostream<_CharT, _Traits>& __os,
1879 const subtract_with_carry_engine<_UI, _W, _S, _R>& __x);
1880
1881 template <class _CharT, class _Traits,
1882 class _UI, size_t _W, size_t _S, size_t _R>
1883 friend
1884 basic_istream<_CharT, _Traits>&
1885 operator>>(basic_istream<_CharT, _Traits>& __is,
1886 subtract_with_carry_engine<_UI, _W, _S, _R>& __x);
1887
1888private:
1889
1890 void seed(result_type __sd, integral_constant<unsigned, 1>);
1891 void seed(result_type __sd, integral_constant<unsigned, 2>);
1892 template<class _Sseq>
1893 void __seed(_Sseq& __q, integral_constant<unsigned, 1>);
1894 template<class _Sseq>
1895 void __seed(_Sseq& __q, integral_constant<unsigned, 2>);
1896};
1897
1898template<class _UIntType, size_t __w, size_t __s, size_t __r>
1899void
1900subtract_with_carry_engine<_UIntType, __w, __s, __r>::seed(result_type __sd,
1901 integral_constant<unsigned, 1>)
1902{
1903 linear_congruential_engine<result_type, 40014u, 0u, 2147483563u>
1904 __e(__sd == 0u ? default_seed : __sd);
1905 for (size_t __i = 0; __i < __r; ++__i)
1906 __x_[__i] = static_cast<result_type>(__e() & _Max);
1907 __c_ = __x_[__r-1] == 0;
1908 __i_ = 0;
1909}
1910
1911template<class _UIntType, size_t __w, size_t __s, size_t __r>
1912void
1913subtract_with_carry_engine<_UIntType, __w, __s, __r>::seed(result_type __sd,
1914 integral_constant<unsigned, 2>)
1915{
1916 linear_congruential_engine<result_type, 40014u, 0u, 2147483563u>
1917 __e(__sd == 0u ? default_seed : __sd);
1918 for (size_t __i = 0; __i < __r; ++__i)
1919 __x_[__i] = static_cast<result_type>(
1920 (__e() + ((uint64_t)__e() << 32)) & _Max);
1921 __c_ = __x_[__r-1] == 0;
1922 __i_ = 0;
1923}
1924
1925template<class _UIntType, size_t __w, size_t __s, size_t __r>
1926template<class _Sseq>
1927void
1928subtract_with_carry_engine<_UIntType, __w, __s, __r>::__seed(_Sseq& __q,
1929 integral_constant<unsigned, 1>)
1930{
1931 const unsigned __k = 1;
1932 uint32_t __ar[__r * __k];
1933 __q.generate(__ar, __ar + __r * __k);
1934 for (size_t __i = 0; __i < __r; ++__i)
1935 __x_[__i] = static_cast<result_type>(__ar[__i] & _Max);
1936 __c_ = __x_[__r-1] == 0;
1937 __i_ = 0;
1938}
1939
1940template<class _UIntType, size_t __w, size_t __s, size_t __r>
1941template<class _Sseq>
1942void
1943subtract_with_carry_engine<_UIntType, __w, __s, __r>::__seed(_Sseq& __q,
1944 integral_constant<unsigned, 2>)
1945{
1946 const unsigned __k = 2;
1947 uint32_t __ar[__r * __k];
1948 __q.generate(__ar, __ar + __r * __k);
1949 for (size_t __i = 0; __i < __r; ++__i)
1950 __x_[__i] = static_cast<result_type>(
1951 (__ar[2 * __i] + ((uint64_t)__ar[2 * __i + 1] << 32)) & _Max);
1952 __c_ = __x_[__r-1] == 0;
1953 __i_ = 0;
1954}
1955
1956template<class _UIntType, size_t __w, size_t __s, size_t __r>
1957_UIntType
1958subtract_with_carry_engine<_UIntType, __w, __s, __r>::operator()()
1959{
1960 const result_type& __xs = __x_[(__i_ + (__r - __s)) % __r];
1961 result_type& __xr = __x_[__i_];
1962 result_type __new_c = __c_ == 0 ? __xs < __xr : __xs != 0 ? __xs <= __xr : 1;
1963 __xr = (__xs - __xr - __c_) & _Max;
1964 __c_ = __new_c;
1965 __i_ = (__i_ + 1) % __r;
1966 return __xr;
1967}
1968
1969template<class _UI, size_t _W, size_t _S, size_t _R>
1970bool
1971operator==(
1972 const subtract_with_carry_engine<_UI, _W, _S, _R>& __x,
1973 const subtract_with_carry_engine<_UI, _W, _S, _R>& __y)
1974{
1975 if (__x.__c_ != __y.__c_)
1976 return false;
1977 if (__x.__i_ == __y.__i_)
1978 return _STD::equal(__x.__x_, __x.__x_ + _R, __y.__x_);
1979 if (__x.__i_ == 0 || __y.__i_ == 0)
1980 {
1981 size_t __j = _STD::min(_R - __x.__i_, _R - __y.__i_);
1982 if (!_STD::equal(__x.__x_ + __x.__i_, __x.__x_ + __x.__i_ + __j,
1983 __y.__x_ + __y.__i_))
1984 return false;
1985 if (__x.__i_ == 0)
1986 return _STD::equal(__x.__x_ + __j, __x.__x_ + _R, __y.__x_);
1987 return _STD::equal(__x.__x_, __x.__x_ + (_R - __j), __y.__x_ + __j);
1988 }
1989 if (__x.__i_ < __y.__i_)
1990 {
1991 size_t __j = _R - __y.__i_;
1992 if (!_STD::equal(__x.__x_ + __x.__i_, __x.__x_ + (__x.__i_ + __j),
1993 __y.__x_ + __y.__i_))
1994 return false;
1995 if (!_STD::equal(__x.__x_ + (__x.__i_ + __j), __x.__x_ + _R,
1996 __y.__x_))
1997 return false;
1998 return _STD::equal(__x.__x_, __x.__x_ + __x.__i_,
1999 __y.__x_ + (_R - (__x.__i_ + __j)));
2000 }
2001 size_t __j = _R - __x.__i_;
2002 if (!_STD::equal(__y.__x_ + __y.__i_, __y.__x_ + (__y.__i_ + __j),
2003 __x.__x_ + __x.__i_))
2004 return false;
2005 if (!_STD::equal(__y.__x_ + (__y.__i_ + __j), __y.__x_ + _R,
2006 __x.__x_))
2007 return false;
2008 return _STD::equal(__y.__x_, __y.__x_ + __y.__i_,
2009 __x.__x_ + (_R - (__y.__i_ + __j)));
2010}
2011
2012template<class _UI, size_t _W, size_t _S, size_t _R>
2013inline
2014bool
2015operator!=(
2016 const subtract_with_carry_engine<_UI, _W, _S, _R>& __x,
2017 const subtract_with_carry_engine<_UI, _W, _S, _R>& __y)
2018{
2019 return !(__x == __y);
2020}
2021
2022template <class _CharT, class _Traits,
2023 class _UI, size_t _W, size_t _S, size_t _R>
2024basic_ostream<_CharT, _Traits>&
2025operator<<(basic_ostream<_CharT, _Traits>& __os,
2026 const subtract_with_carry_engine<_UI, _W, _S, _R>& __x)
2027{
2028 __save_flags<_CharT, _Traits> _(__os);
2029 __os.flags(ios_base::dec | ios_base::left);
2030 _CharT __sp = __os.widen(' ');
2031 __os.fill(__sp);
2032 __os << __x.__x_[__x.__i_];
2033 for (size_t __j = __x.__i_ + 1; __j < _R; ++__j)
2034 __os << __sp << __x.__x_[__j];
2035 for (size_t __j = 0; __j < __x.__i_; ++__j)
2036 __os << __sp << __x.__x_[__j];
2037 __os << __sp << __x.__c_;
2038 return __os;
2039}
2040
2041template <class _CharT, class _Traits,
2042 class _UI, size_t _W, size_t _S, size_t _R>
2043basic_istream<_CharT, _Traits>&
2044operator>>(basic_istream<_CharT, _Traits>& __is,
2045 subtract_with_carry_engine<_UI, _W, _S, _R>& __x)
2046{
2047 __save_flags<_CharT, _Traits> _(__is);
2048 __is.flags(ios_base::dec | ios_base::skipws);
2049 _UI __t[_R+1];
2050 for (size_t __i = 0; __i < _R+1; ++__i)
2051 __is >> __t[__i];
2052 if (!__is.fail())
2053 {
2054 for (size_t __i = 0; __i < _R; ++__i)
2055 __x.__x_[__i] = __t[__i];
2056 __x.__c_ = __t[_R];
2057 __x.__i_ = 0;
2058 }
2059 return __is;
2060}
2061
2062typedef subtract_with_carry_engine<uint_fast32_t, 24, 10, 24> ranlux24_base;
2063typedef subtract_with_carry_engine<uint_fast64_t, 48, 5, 12> ranlux48_base;
2064
2065// discard_block_engine
2066
2067template<class _Engine, size_t __p, size_t __r>
2068class discard_block_engine
2069{
2070 _Engine __e_;
2071 int __n_;
2072
2073 static_assert( 0 < __r, "discard_block_engine invalid parameters");
2074 static_assert(__r <= __p, "discard_block_engine invalid parameters");
2075public:
2076 // types
2077 typedef typename _Engine::result_type result_type;
2078
2079 // engine characteristics
2080 static const/*expr*/ size_t block_size = __p;
2081 static const/*expr*/ size_t used_block = __r;
2082
2083 // Temporary work around for lack of constexpr
2084 static const result_type _Min = _Engine::_Min;
2085 static const result_type _Max = _Engine::_Max;
2086
2087 static const/*expr*/ result_type min() { return _Engine::min(); }
2088 static const/*expr*/ result_type max() { return _Engine::max(); }
2089
2090 // constructors and seeding functions
2091 discard_block_engine() : __n_(0) {}
2092// explicit discard_block_engine(const _Engine& __e);
2093// explicit discard_block_engine(_Engine&& __e);
2094 explicit discard_block_engine(result_type __sd) : __e_(__sd), __n_(0) {}
2095 template<class _Sseq> explicit discard_block_engine(_Sseq& __q)
2096 : __e_(__q), __n_(0) {}
2097 void seed() {__e_.seed(); __n_ = 0;}
2098 void seed(result_type __sd) {__e_.seed(__sd); __n_ = 0;}
2099 template<class _Sseq> void seed(_Sseq& __q) {__e_.seed(__q); __n_ = 0;}
2100
2101 // generating functions
2102 result_type operator()();
2103 void discard(unsigned long long __z) {for (; __z; --__z) operator()();}
2104
2105 // property functions
2106 const _Engine& base() const {return __e_;}
2107
2108 template<class _Eng, size_t _P, size_t _R>
2109 friend
2110 bool
2111 operator==(
2112 const discard_block_engine<_Eng, _P, _R>& __x,
2113 const discard_block_engine<_Eng, _P, _R>& __y);
2114
2115 template<class _Eng, size_t _P, size_t _R>
2116 friend
2117 bool
2118 operator!=(
2119 const discard_block_engine<_Eng, _P, _R>& __x,
2120 const discard_block_engine<_Eng, _P, _R>& __y);
2121
2122 template <class _CharT, class _Traits,
2123 class _Eng, size_t _P, size_t _R>
2124 friend
2125 basic_ostream<_CharT, _Traits>&
2126 operator<<(basic_ostream<_CharT, _Traits>& __os,
2127 const discard_block_engine<_Eng, _P, _R>& __x);
2128
2129 template <class _CharT, class _Traits,
2130 class _Eng, size_t _P, size_t _R>
2131 friend
2132 basic_istream<_CharT, _Traits>&
2133 operator>>(basic_istream<_CharT, _Traits>& __is,
2134 discard_block_engine<_Eng, _P, _R>& __x);
2135};
2136
2137template<class _Engine, size_t __p, size_t __r>
2138typename discard_block_engine<_Engine, __p, __r>::result_type
2139discard_block_engine<_Engine, __p, __r>::operator()()
2140{
2141 if (__n_ >= __r)
2142 {
2143 __e_.discard(__p - __r);
2144 __n_ = 0;
2145 }
2146 ++__n_;
2147 return __e_();
2148}
2149
2150template<class _Eng, size_t _P, size_t _R>
2151inline
2152bool
2153operator==(const discard_block_engine<_Eng, _P, _R>& __x,
2154 const discard_block_engine<_Eng, _P, _R>& __y)
2155{
2156 return __x.__n_ == __y.__n_ && __x.__e_ == __y.__e_;
2157}
2158
2159template<class _Eng, size_t _P, size_t _R>
2160inline
2161bool
2162operator!=(const discard_block_engine<_Eng, _P, _R>& __x,
2163 const discard_block_engine<_Eng, _P, _R>& __y)
2164{
2165 return !(__x == __y);
2166}
2167
2168template <class _CharT, class _Traits,
2169 class _Eng, size_t _P, size_t _R>
2170basic_ostream<_CharT, _Traits>&
2171operator<<(basic_ostream<_CharT, _Traits>& __os,
2172 const discard_block_engine<_Eng, _P, _R>& __x)
2173{
2174 __save_flags<_CharT, _Traits> _(__os);
2175 __os.flags(ios_base::dec | ios_base::left);
2176 _CharT __sp = __os.widen(' ');
2177 __os.fill(__sp);
2178 return __os << __x.__e_ << __sp << __x.__n_;
2179}
2180
2181template <class _CharT, class _Traits,
2182 class _Eng, size_t _P, size_t _R>
2183basic_istream<_CharT, _Traits>&
2184operator>>(basic_istream<_CharT, _Traits>& __is,
2185 discard_block_engine<_Eng, _P, _R>& __x)
2186{
2187 __save_flags<_CharT, _Traits> _(__is);
2188 __is.flags(ios_base::dec | ios_base::skipws);
2189 _Eng __e;
2190 int __n;
2191 __is >> __e >> __n;
2192 if (!__is.fail())
2193 {
2194 __x.__e_ = __e;
2195 __x.__n_ = __n;
2196 }
2197 return __is;
2198}
2199
2200typedef discard_block_engine<ranlux24_base, 223, 23> ranlux24;
2201typedef discard_block_engine<ranlux48_base, 389, 11> ranlux48;
2202
2203// independent_bits_engine
2204
2205template <unsigned long long _X, size_t _R>
2206struct __log2_imp
2207{
2208 static const size_t value = _X & ((unsigned long long)(1) << _R) ? _R
2209 : __log2_imp<_X, _R - 1>::value;
2210};
2211
2212template <unsigned long long _X>
2213struct __log2_imp<_X, 0>
2214{
2215 static const size_t value = 0;
2216};
2217
2218template <size_t _R>
2219struct __log2_imp<0, _R>
2220{
2221 static const size_t value = _R + 1;
2222};
2223
2224template <class _UI, _UI _X>
2225struct __log2
2226{
2227 static const size_t value = __log2_imp<_X,
2228 sizeof(_UI) * __CHAR_BIT__ - 1>::value;
2229};
2230
2231template<class _Engine, size_t __w, class _UIntType>
2232class independent_bits_engine
2233{
2234 template <class _UI, _UI _R0, size_t _W, size_t _M>
2235 class __get_n
2236 {
2237 static const size_t _Dt = numeric_limits<_UI>::digits;
2238 static const size_t _N = _W / _M + (_W % _M != 0);
2239 static const size_t _W0 = _W / _N;
2240 static const _UI _Y0 = _W0 >= _Dt ? 0 : (_R0 >> _W0) << _W0;
2241 public:
2242 static const size_t value = _R0 - _Y0 > _Y0 / _N ? _N + 1 : _N;
2243 };
2244public:
2245 // types
2246 typedef _UIntType result_type;
2247
2248private:
2249 _Engine __e_;
2250
2251 static const result_type _Dt = numeric_limits<result_type>::digits;
2252 static_assert( 0 < __w, "independent_bits_engine invalid parameters");
2253 static_assert(__w <= _Dt, "independent_bits_engine invalid parameters");
2254
2255 typedef typename _Engine::result_type _Engine_result_type;
2256 typedef typename conditional
2257 <
2258 sizeof(_Engine_result_type) <= sizeof(result_type),
2259 result_type,
2260 _Engine_result_type
2261 >::type _Working_result_type;
2262 // Temporary work around for lack of constexpr
2263 static const _Working_result_type _R = _Engine::_Max - _Engine::_Min
2264 + _Working_result_type(1);
2265 static const size_t __m = __log2<_Working_result_type, _R>::value;
2266 static const size_t __n = __get_n<_Working_result_type, _R, __w, __m>::value;
2267 static const size_t __w0 = __w / __n;
2268 static const size_t __n0 = __n - __w % __n;
2269 static const size_t _WDt = numeric_limits<_Working_result_type>::digits;
2270 static const size_t _EDt = numeric_limits<_Engine_result_type>::digits;
2271 static const _Working_result_type __y0 = __w0 >= _WDt ? 0 :
2272 (_R >> __w0) << __w0;
2273 static const _Working_result_type __y1 = __w0 >= _WDt - 1 ? 0 :
2274 (_R >> (__w0+1)) << (__w0+1);
2275 static const _Engine_result_type __mask0 = __w0 > 0 ?
2276 _Engine_result_type(~0) >> (_EDt - __w0) :
2277 _Engine_result_type(0);
2278 static const _Engine_result_type __mask1 = __w0 < _EDt - 1 ?
2279 _Engine_result_type(~0) >> (_EDt - (__w0 + 1)) :
2280 _Engine_result_type(~0);
2281public:
2282 static const result_type _Min = 0;
2283 static const result_type _Max = __w == _Dt ? result_type(~0) :
2284 (result_type(1) << __w) - result_type(1);
2285 static_assert(_Min < _Max, "independent_bits_engine invalid parameters");
2286
2287 // engine characteristics
2288 static const/*expr*/ result_type min() { return _Min; }
2289 static const/*expr*/ result_type max() { return _Max; }
2290
2291 // constructors and seeding functions
2292 independent_bits_engine() {}
2293// explicit independent_bits_engine(const _Engine& __e);
2294// explicit independent_bits_engine(_Engine&& __e);
2295 explicit independent_bits_engine(result_type __sd) : __e_(__sd) {}
2296 template<class _Sseq> explicit independent_bits_engine(_Sseq& __q)
2297 : __e_(__q) {}
2298 void seed() {__e_.seed();}
2299 void seed(result_type __sd) {__e_.seed(__sd);}
2300 template<class _Sseq> void seed(_Sseq& __q) {__e_.seed(__q);}
2301
2302 // generating functions
2303 result_type operator()() {return __eval(integral_constant<bool, _R != 0>());}
2304 void discard(unsigned long long __z) {for (; __z; --__z) operator()();}
2305
2306 // property functions
2307 const _Engine& base() const {return __e_;}
2308
2309 template<class _Eng, size_t _W, class _UI>
2310 friend
2311 bool
2312 operator==(
2313 const independent_bits_engine<_Eng, _W, _UI>& __x,
2314 const independent_bits_engine<_Eng, _W, _UI>& __y);
2315
2316 template<class _Eng, size_t _W, class _UI>
2317 friend
2318 bool
2319 operator!=(
2320 const independent_bits_engine<_Eng, _W, _UI>& __x,
2321 const independent_bits_engine<_Eng, _W, _UI>& __y);
2322
2323 template <class _CharT, class _Traits,
2324 class _Eng, size_t _W, class _UI>
2325 friend
2326 basic_ostream<_CharT, _Traits>&
2327 operator<<(basic_ostream<_CharT, _Traits>& __os,
2328 const independent_bits_engine<_Eng, _W, _UI>& __x);
2329
2330 template <class _CharT, class _Traits,
2331 class _Eng, size_t _W, class _UI>
2332 friend
2333 basic_istream<_CharT, _Traits>&
2334 operator>>(basic_istream<_CharT, _Traits>& __is,
2335 independent_bits_engine<_Eng, _W, _UI>& __x);
2336
2337private:
2338 result_type __eval(false_type);
2339 result_type __eval(true_type);
2340
2341 template <size_t __count>
2342 static
2343 typename enable_if
2344 <
2345 __count < _Dt,
2346 result_type
2347 >::type
2348 __lshift(result_type __x) {return __x << __count;}
2349
2350 template <size_t __count>
2351 static
2352 typename enable_if
2353 <
2354 (__count >= _Dt),
2355 result_type
2356 >::type
2357 __lshift(result_type __x) {return result_type(0);}
2358};
2359
2360template<class _Engine, size_t __w, class _UIntType>
2361inline
2362_UIntType
2363independent_bits_engine<_Engine, __w, _UIntType>::__eval(false_type)
2364{
2365 return static_cast<result_type>(__e_() & __mask0);
2366}
2367
2368template<class _Engine, size_t __w, class _UIntType>
2369_UIntType
2370independent_bits_engine<_Engine, __w, _UIntType>::__eval(true_type)
2371{
2372 result_type _S = 0;
2373 for (size_t __k = 0; __k < __n0; ++__k)
2374 {
2375 _Engine_result_type __u;
2376 do
2377 {
2378 __u = __e_() - _Engine::min();
2379 } while (__u >= __y0);
2380 _S = static_cast<result_type>(__lshift<__w0>(_S) + (__u & __mask0));
2381 }
2382 for (size_t __k = __n0; __k < __n; ++__k)
2383 {
2384 _Engine_result_type __u;
2385 do
2386 {
2387 __u = __e_() - _Engine::min();
2388 } while (__u >= __y1);
2389 _S = static_cast<result_type>(__lshift<__w0+1>(_S) + (__u & __mask1));
2390 }
2391 return _S;
2392}
2393
2394template<class _Eng, size_t _W, class _UI>
2395inline
2396bool
2397operator==(
2398 const independent_bits_engine<_Eng, _W, _UI>& __x,
2399 const independent_bits_engine<_Eng, _W, _UI>& __y)
2400{
2401 return __x.base() == __y.base();
2402}
2403
2404template<class _Eng, size_t _W, class _UI>
2405inline
2406bool
2407operator!=(
2408 const independent_bits_engine<_Eng, _W, _UI>& __x,
2409 const independent_bits_engine<_Eng, _W, _UI>& __y)
2410{
2411 return !(__x == __y);
2412}
2413
2414template <class _CharT, class _Traits,
2415 class _Eng, size_t _W, class _UI>
2416basic_ostream<_CharT, _Traits>&
2417operator<<(basic_ostream<_CharT, _Traits>& __os,
2418 const independent_bits_engine<_Eng, _W, _UI>& __x)
2419{
2420 return __os << __x.base();
2421}
2422
2423template <class _CharT, class _Traits,
2424 class _Eng, size_t _W, class _UI>
2425basic_istream<_CharT, _Traits>&
2426operator>>(basic_istream<_CharT, _Traits>& __is,
2427 independent_bits_engine<_Eng, _W, _UI>& __x)
2428{
2429 _Eng __e;
2430 __is >> __e;
2431 if (!__is.fail())
2432 __x.__e_ = __e;
2433 return __is;
2434}
2435
2436// shuffle_order_engine
2437
2438template <uint64_t _Xp, uint64_t _Yp>
2439struct __ugcd
2440{
2441 static const uint64_t value = __ugcd<_Yp, _Xp % _Yp>::value;
2442};
2443
2444template <uint64_t _Xp>
2445struct __ugcd<_Xp, 0>
2446{
2447 static const uint64_t value = _Xp;
2448};
2449
2450template <uint64_t _N, uint64_t _D>
2451class __uratio
2452{
2453 static_assert(_D != 0, "__uratio divide by 0");
2454 static const uint64_t __gcd = __ugcd<_N, _D>::value;
2455public:
2456 static const uint64_t num = _N / __gcd;
2457 static const uint64_t den = _D / __gcd;
2458
2459 typedef __uratio<num, den> type;
2460};
2461
2462template<class _Engine, size_t __k>
2463class shuffle_order_engine
2464{
2465 static_assert(0 < __k, "shuffle_order_engine invalid parameters");
2466public:
2467 // types
2468 typedef typename _Engine::result_type result_type;
2469
2470private:
2471 _Engine __e_;
2472 result_type _V_[__k];
2473 result_type _Y_;
2474
2475public:
2476 // engine characteristics
2477 static const/*expr*/ size_t table_size = __k;
2478
2479 static const result_type _Min = _Engine::_Min;
2480 static const result_type _Max = _Engine::_Max;
2481 static_assert(_Min < _Max, "shuffle_order_engine invalid parameters");
2482 static const/*expr*/ result_type min() { return _Min; }
2483 static const/*expr*/ result_type max() { return _Max; }
2484
2485 static const unsigned long long _R = _Max - _Min + 1ull;
2486
2487 // constructors and seeding functions
2488 shuffle_order_engine() {__init();}
2489// explicit shuffle_order_engine(const _Engine& __e);
2490// explicit shuffle_order_engine(_Engine&& e);
2491 explicit shuffle_order_engine(result_type __sd) : __e_(__sd) {__init();}
2492 template<class _Sseq> explicit shuffle_order_engine(_Sseq& __q)
2493 : __e_(__q) {__init();}
2494 void seed() {__e_.seed(); __init();}
2495 void seed(result_type __sd) {__e_.seed(__sd); __init();}
2496 template<class _Sseq> void seed(_Sseq& __q) {__e_.seed(__q); __init();}
2497
2498 // generating functions
2499 result_type operator()() {return __eval(integral_constant<bool, _R != 0>());}
2500 void discard(unsigned long long __z) {for (; __z; --__z) operator()();}
2501
2502 // property functions
2503 const _Engine& base() const {return __e_;}
2504
2505private:
2506 template<class _Eng, size_t _K>
2507 friend
2508 bool
2509 operator==(
2510 const shuffle_order_engine<_Eng, _K>& __x,
2511 const shuffle_order_engine<_Eng, _K>& __y);
2512
2513 template<class _Eng, size_t _K>
2514 friend
2515 bool
2516 operator!=(
2517 const shuffle_order_engine<_Eng, _K>& __x,
2518 const shuffle_order_engine<_Eng, _K>& __y);
2519
2520 template <class _CharT, class _Traits,
2521 class _Eng, size_t _K>
2522 friend
2523 basic_ostream<_CharT, _Traits>&
2524 operator<<(basic_ostream<_CharT, _Traits>& __os,
2525 const shuffle_order_engine<_Eng, _K>& __x);
2526
2527 template <class _CharT, class _Traits,
2528 class _Eng, size_t _K>
2529 friend
2530 basic_istream<_CharT, _Traits>&
2531 operator>>(basic_istream<_CharT, _Traits>& __is,
2532 shuffle_order_engine<_Eng, _K>& __x);
2533
2534 void __init()
2535 {
2536 for (size_t __i = 0; __i < __k; ++__i)
2537 _V_[__i] = __e_();
2538 _Y_ = __e_();
2539 }
2540
2541 result_type __eval(false_type) {return __eval2(integral_constant<bool, __k & 1>());}
2542 result_type __eval(true_type) {return __eval(__uratio<__k, _R>());}
2543
2544 result_type __eval2(false_type) {return __eval(__uratio<__k/2, 0x8000000000000000ull>());}
2545 result_type __eval2(true_type) {return __evalf<__k, 0>();}
2546
2547 template <uint64_t _N, uint64_t _D>
2548 typename enable_if
2549 <
2550 (__uratio<_N, _D>::num > 0xFFFFFFFFFFFFFFFFull / (_Max - _Min)),
2551 result_type
2552 >::type
2553 __eval(__uratio<_N, _D>)
2554 {return __evalf<__uratio<_N, _D>::num, __uratio<_N, _D>::den>();}
2555
2556 template <uint64_t _N, uint64_t _D>
2557 typename enable_if
2558 <
2559 __uratio<_N, _D>::num <= 0xFFFFFFFFFFFFFFFFull / (_Max - _Min),
2560 result_type
2561 >::type
2562 __eval(__uratio<_N, _D>)
2563 {
2564 const size_t __j = static_cast<size_t>(__uratio<_N, _D>::num * (_Y_ - _Min)
2565 / __uratio<_N, _D>::den);
2566 _Y_ = _V_[__j];
2567 _V_[__j] = __e_();
2568 return _Y_;
2569 }
2570
2571 template <uint64_t __n, uint64_t __d>
2572 result_type __evalf()
2573 {
2574 const double _F = __d == 0 ?
2575 __n / (2. * 0x8000000000000000ull) :
2576 __n / (double)__d;
2577 const size_t __j = static_cast<size_t>(_F * (_Y_ - _Min));
2578 _Y_ = _V_[__j];
2579 _V_[__j] = __e_();
2580 return _Y_;
2581 }
2582};
2583
2584template<class _Eng, size_t _K>
2585bool
2586operator==(
2587 const shuffle_order_engine<_Eng, _K>& __x,
2588 const shuffle_order_engine<_Eng, _K>& __y)
2589{
2590 return __x._Y_ == __y._Y_ && _STD::equal(__x._V_, __x._V_ + _K, __y._V_) &&
2591 __x.__e_ == __y.__e_;
2592}
2593
2594template<class _Eng, size_t _K>
2595inline
2596bool
2597operator!=(
2598 const shuffle_order_engine<_Eng, _K>& __x,
2599 const shuffle_order_engine<_Eng, _K>& __y)
2600{
2601 return !(__x == __y);
2602}
2603
2604template <class _CharT, class _Traits,
2605 class _Eng, size_t _K>
2606basic_ostream<_CharT, _Traits>&
2607operator<<(basic_ostream<_CharT, _Traits>& __os,
2608 const shuffle_order_engine<_Eng, _K>& __x)
2609{
2610 __save_flags<_CharT, _Traits> _(__os);
2611 __os.flags(ios_base::dec | ios_base::left);
2612 _CharT __sp = __os.widen(' ');
2613 __os.fill(__sp);
2614 __os << __x.__e_ << __sp << __x._V_[0];
2615 for (size_t __i = 1; __i < _K; ++__i)
2616 __os << __sp << __x._V_[__i];
2617 return __os << __sp << __x._Y_;
2618}
2619
2620template <class _CharT, class _Traits,
2621 class _Eng, size_t _K>
2622basic_istream<_CharT, _Traits>&
2623operator>>(basic_istream<_CharT, _Traits>& __is,
2624 shuffle_order_engine<_Eng, _K>& __x)
2625{
2626 typedef typename shuffle_order_engine<_Eng, _K>::result_type result_type;
2627 __save_flags<_CharT, _Traits> _(__is);
2628 __is.flags(ios_base::dec | ios_base::skipws);
2629 _Eng __e;
2630 result_type _V[_K+1];
2631 __is >> __e;
2632 for (size_t __i = 0; __i < _K+1; ++__i)
2633 __is >> _V[__i];
2634 if (!__is.fail())
2635 {
2636 __x.__e_ = __e;
2637 for (size_t __i = 0; __i < _K; ++__i)
2638 __x._V_[__i] = _V[__i];
2639 __x._Y_ = _V[_K];
2640 }
2641 return __is;
2642}
2643
2644typedef shuffle_order_engine<minstd_rand0, 256> knuth_b;
2645
2646// random_device
2647
2648class random_device
2649{
2650 int __f_;
2651public:
2652 // types
2653 typedef unsigned result_type;
2654
2655 // generator characteristics
2656 static const result_type _Min = 0;
2657 static const result_type _Max = 0xFFFFFFFFu;
2658
2659 static const/*expr*/ result_type min() { return _Min;}
2660 static const/*expr*/ result_type max() { return _Max;}
2661
2662 // constructors
2663 explicit random_device(const string& __token = "/dev/urandom");
2664 ~random_device();
2665
2666 // generating functions
2667 result_type operator()();
2668
2669 // property functions
2670 double entropy() const;
2671
2672private:
2673 // no copy functions
2674 random_device(const random_device&); // = delete;
2675 random_device& operator=(const random_device&); // = delete;
2676};
2677
2678// seed_seq
2679
2680class seed_seq
2681{
2682public:
2683 // types
2684 typedef uint32_t result_type;
2685
2686private:
2687 vector<result_type> __v_;
2688
2689 template<class _InputIterator>
2690 void init(_InputIterator __first, _InputIterator __last);
2691public:
2692 // constructors
2693 seed_seq() {}
2694 template<class _Tp>
2695 seed_seq(initializer_list<_Tp> __il) {init(__il.begin(), __il.end());}
2696
2697 template<class _InputIterator>
2698 seed_seq(_InputIterator __first, _InputIterator __last)
2699 {init(__first, __last);}
2700
2701 // generating functions
2702 template<class _RandomAccessIterator>
2703 void generate(_RandomAccessIterator __first, _RandomAccessIterator __last);
2704
2705 // property functions
2706 size_t size() const {return __v_.size();}
2707 template<class _OutputIterator>
2708 void param(_OutputIterator __dest) const
2709 {_STD::copy(__v_.begin(), __v_.end(), __dest);}
2710
2711private:
2712 // no copy functions
2713 seed_seq(const seed_seq&); // = delete;
2714 void operator=(const seed_seq&); // = delete;
2715
2716 static result_type _T(result_type __x) {return __x ^ (__x >> 27);}
2717};
2718
2719template<class _InputIterator>
2720void
2721seed_seq::init(_InputIterator __first, _InputIterator __last)
2722{
2723 for (_InputIterator __s = __first; __s != __last; ++__s)
2724 __v_.push_back(*__s & 0xFFFFFFFF);
2725}
2726
2727template<class _RandomAccessIterator>
2728void
2729seed_seq::generate(_RandomAccessIterator __first, _RandomAccessIterator __last)
2730{
2731 if (__first != __last)
2732 {
2733 _STD::fill(__first, __last, 0x8b8b8b8b);
2734 const size_t __n = static_cast<size_t>(__last - __first);
2735 const size_t __s = __v_.size();
2736 const size_t __t = (__n >= 623) ? 11
2737 : (__n >= 68) ? 7
2738 : (__n >= 39) ? 5
2739 : (__n >= 7) ? 3
2740 : (__n - 1) / 2;
2741 const size_t __p = (__n - __t) / 2;
2742 const size_t __q = __p + __t;
2743 const size_t __m = _STD::max(__s + 1, __n);
2744 // __k = 0;
2745 {
2746 result_type __r = 1664525 * _T(__first[0] ^ __first[__p]
2747 ^ __first[__n - 1]);
2748 __first[__p] += __r;
2749 __r += __s;
2750 __first[__q] += __r;
2751 __first[0] = __r;
2752 }
2753 for (size_t __k = 1; __k <= __s; ++__k)
2754 {
2755 const size_t __kmodn = __k % __n;
2756 const size_t __kpmodn = (__k + __p) % __n;
2757 result_type __r = 1664525 * _T(__first[__kmodn] ^ __first[__kpmodn]
2758 ^ __first[(__k - 1) % __n]);
2759 __first[__kpmodn] += __r;
2760 __r += __kmodn + __v_[__k-1];
2761 __first[(__k + __q) % __n] += __r;
2762 __first[__kmodn] = __r;
2763 }
2764 for (size_t __k = __s + 1; __k < __m; ++__k)
2765 {
2766 const size_t __kmodn = __k % __n;
2767 const size_t __kpmodn = (__k + __p) % __n;
2768 result_type __r = 1664525 * _T(__first[__kmodn] ^ __first[__kpmodn]
2769 ^ __first[(__k - 1) % __n]);
2770 __first[__kpmodn] += __r;
2771 __r += __kmodn;
2772 __first[(__k + __q) % __n] += __r;
2773 __first[__kmodn] = __r;
2774 }
2775 for (size_t __k = __m; __k < __m + __n; ++__k)
2776 {
2777 const size_t __kmodn = __k % __n;
2778 const size_t __kpmodn = (__k + __p) % __n;
2779 result_type __r = 1566083941 * _T(__first[__kmodn] +
2780 __first[__kpmodn] +
2781 __first[(__k - 1) % __n]);
2782 __first[__kpmodn] ^= __r;
2783 __r -= __kmodn;
2784 __first[(__k + __q) % __n] ^= __r;
2785 __first[__kmodn] = __r;
2786 }
2787 }
2788}
2789
Howard Hinnant30a840f2010-05-12 17:08:57 +00002790// generate_canonical
2791
Howard Hinnantbc8d3f92010-05-11 19:42:16 +00002792template<class _RealType, size_t __bits, class _URNG>
2793_RealType
2794generate_canonical(_URNG& __g)
2795{
2796 const size_t _Dt = numeric_limits<_RealType>::digits;
2797 const size_t __b = _Dt < __bits ? _Dt : __bits;
2798 const size_t __logR = __log2<uint64_t, _URNG::_Max - _URNG::_Min + uint64_t(1)>::value;
2799 const size_t __k = __b / __logR + (__b % __logR != 0) + (__b == 0);
2800 const _RealType _R = _URNG::_Max - _URNG::_Min + _RealType(1);
2801 _RealType __base = _R;
2802 _RealType _S = __g() - _URNG::_Min;
2803 for (size_t __i = 1; __i < __k; ++__i, __base *= _R)
2804 _S += (__g() - _URNG::_Min) * __base;
2805 return _S / __base;
2806}
2807
2808// __independent_bits_engine
2809
2810template<class _Engine, class _UIntType>
2811class __independent_bits_engine
2812{
2813public:
2814 // types
2815 typedef _UIntType result_type;
2816
2817private:
2818 typedef typename _Engine::result_type _Engine_result_type;
2819 typedef typename conditional
2820 <
2821 sizeof(_Engine_result_type) <= sizeof(result_type),
2822 result_type,
2823 _Engine_result_type
2824 >::type _Working_result_type;
2825
2826 _Engine& __e_;
2827 size_t __w_;
2828 size_t __w0_;
2829 size_t __n_;
2830 size_t __n0_;
2831 _Working_result_type __y0_;
2832 _Working_result_type __y1_;
2833 _Engine_result_type __mask0_;
2834 _Engine_result_type __mask1_;
2835
2836 static const _Working_result_type _R = _Engine::_Max - _Engine::_Min
2837 + _Working_result_type(1);
2838 static const size_t __m = __log2<_Working_result_type, _R>::value;
2839 static const size_t _WDt = numeric_limits<_Working_result_type>::digits;
2840 static const size_t _EDt = numeric_limits<_Engine_result_type>::digits;
2841
2842public:
2843 // constructors and seeding functions
2844 __independent_bits_engine(_Engine& __e, size_t __w);
2845
2846 // generating functions
2847 result_type operator()() {return __eval(integral_constant<bool, _R != 0>());}
2848
2849private:
2850 result_type __eval(false_type);
2851 result_type __eval(true_type);
2852};
2853
2854template<class _Engine, class _UIntType>
2855__independent_bits_engine<_Engine, _UIntType>
2856 ::__independent_bits_engine(_Engine& __e, size_t __w)
2857 : __e_(__e),
2858 __w_(__w)
2859{
2860 __n_ = __w_ / __m + (__w_ % __m != 0);
2861 __w0_ = __w_ / __n_;
2862 if (_R == 0)
2863 __y0_ = _R;
2864 else if (__w0_ < _WDt)
2865 __y0_ = (_R >> __w0_) << __w0_;
2866 else
2867 __y0_ = 0;
2868 if (_R - __y0_ > __y0_ / __n_)
2869 {
2870 ++__n_;
2871 __w0_ = __w_ / __n_;
2872 if (__w0_ < _WDt)
2873 __y0_ = (_R >> __w0_) << __w0_;
2874 else
2875 __y0_ = 0;
2876 }
2877 __n0_ = __n_ - __w_ % __n_;
2878 if (__w0_ < _WDt - 1)
2879 __y1_ = (_R >> (__w0_ + 1)) << (__w0_ + 1);
2880 else
2881 __y1_ = 0;
2882 __mask0_ = __w0_ > 0 ? _Engine_result_type(~0) >> (_EDt - __w0_) :
2883 _Engine_result_type(0);
2884 __mask1_ = __w0_ < _EDt - 1 ?
2885 _Engine_result_type(~0) >> (_EDt - (__w0_ + 1)) :
2886 _Engine_result_type(~0);
2887}
2888
2889template<class _Engine, class _UIntType>
2890inline
2891_UIntType
2892__independent_bits_engine<_Engine, _UIntType>::__eval(false_type)
2893{
2894 return static_cast<result_type>(__e_() & __mask0_);
2895}
2896
2897template<class _Engine, class _UIntType>
2898_UIntType
2899__independent_bits_engine<_Engine, _UIntType>::__eval(true_type)
2900{
2901 result_type _S = 0;
2902 for (size_t __k = 0; __k < __n0_; ++__k)
2903 {
2904 _Engine_result_type __u;
2905 do
2906 {
2907 __u = __e_() - _Engine::min();
2908 } while (__u >= __y0_);
2909 if (__w0_ < _EDt)
2910 _S <<= __w0_;
2911 else
2912 _S = 0;
2913 _S += __u & __mask0_;
2914 }
2915 for (size_t __k = __n0_; __k < __n_; ++__k)
2916 {
2917 _Engine_result_type __u;
2918 do
2919 {
2920 __u = __e_() - _Engine::min();
2921 } while (__u >= __y1_);
2922 if (__w0_ < _EDt - 1)
2923 _S <<= __w0_ + 1;
2924 else
2925 _S = 0;
2926 _S += __u & __mask1_;
2927 }
2928 return _S;
2929}
2930
Howard Hinnant30a840f2010-05-12 17:08:57 +00002931// uniform_int_distribution
2932
Howard Hinnantbc8d3f92010-05-11 19:42:16 +00002933template<class _IntType = int>
2934class uniform_int_distribution
2935{
2936public:
2937 // types
2938 typedef _IntType result_type;
2939
2940 class param_type
2941 {
2942 result_type __a_;
2943 result_type __b_;
2944 public:
2945 typedef uniform_int_distribution distribution_type;
2946
2947 explicit param_type(result_type __a = 0,
2948 result_type __b = numeric_limits<result_type>::max())
2949 : __a_(__a), __b_(__b) {}
2950
2951 result_type a() const {return __a_;}
2952 result_type b() const {return __b_;}
2953
2954 friend bool operator==(const param_type& __x, const param_type& __y)
2955 {return __x.__a_ == __y.__a_ && __x.__b_ == __y.__b_;}
2956 friend bool operator!=(const param_type& __x, const param_type& __y)
2957 {return !(__x == __y);}
2958 };
2959
2960private:
2961 param_type __p_;
2962
2963public:
2964 // constructors and reset functions
2965 explicit uniform_int_distribution(result_type __a = 0,
2966 result_type __b = numeric_limits<result_type>::max())
2967 : __p_(param_type(__a, __b)) {}
2968 explicit uniform_int_distribution(const param_type& __p) : __p_(__p) {}
2969 void reset() {}
2970
2971 // generating functions
2972 template<class _URNG> result_type operator()(_URNG& __g)
2973 {return (*this)(__g, __p_);}
2974 template<class _URNG> result_type operator()(_URNG& __g, const param_type& __p);
2975
2976 // property functions
2977 result_type a() const {return __p_.a();}
2978 result_type b() const {return __p_.b();}
2979
2980 param_type param() const {return __p_;}
2981 void param(const param_type& __p) {__p_ = __p;}
2982
2983 result_type min() const {return a();}
2984 result_type max() const {return b();}
2985
2986 friend bool operator==(const uniform_int_distribution& __x,
2987 const uniform_int_distribution& __y)
2988 {return __x.__p_ == __y.__p_;}
2989 friend bool operator!=(const uniform_int_distribution& __x,
2990 const uniform_int_distribution& __y)
2991 {return !(__x == __y);}
2992};
2993
2994template<class _IntType>
2995template<class _URNG>
2996typename uniform_int_distribution<_IntType>::result_type
2997uniform_int_distribution<_IntType>::operator()(_URNG& __g, const param_type& __p)
2998{
2999 typedef typename conditional<sizeof(result_type) <= sizeof(uint32_t),
3000 uint32_t, uint64_t>::type _UIntType;
3001 const _UIntType _R = __p.b() - __p.a() + _UIntType(1);
3002 if (_R == 1)
3003 return __p.a();
3004 const size_t _Dt = numeric_limits<_UIntType>::digits;
3005 typedef __independent_bits_engine<_URNG, _UIntType> _Eng;
3006 if (_R == 0)
3007 return static_cast<result_type>(_Eng(__g, _Dt)());
3008 size_t __w = _Dt - __clz(_R) - 1;
3009 if ((_R & (_UIntType(~0) >> (_Dt - __w))) != 0)
3010 ++__w;
3011 _Eng __e(__g, __w);
3012 _UIntType __u;
3013 do
3014 {
3015 __u = __e();
3016 } while (__u >= _R);
3017 return static_cast<result_type>(__u + __p.a());
3018}
3019
3020template <class _CharT, class _Traits, class _IT>
3021basic_ostream<_CharT, _Traits>&
3022operator<<(basic_ostream<_CharT, _Traits>& __os,
3023 const uniform_int_distribution<_IT>& __x)
3024{
3025 __save_flags<_CharT, _Traits> _(__os);
3026 __os.flags(ios_base::dec | ios_base::left);
3027 _CharT __sp = __os.widen(' ');
3028 __os.fill(__sp);
3029 return __os << __x.a() << __sp << __x.b();
3030}
3031
3032template <class _CharT, class _Traits, class _IT>
3033basic_istream<_CharT, _Traits>&
3034operator>>(basic_istream<_CharT, _Traits>& __is,
3035 uniform_int_distribution<_IT>& __x)
3036{
3037 typedef uniform_int_distribution<_IT> _Eng;
3038 typedef typename _Eng::result_type result_type;
3039 typedef typename _Eng::param_type param_type;
3040 __save_flags<_CharT, _Traits> _(__is);
3041 __is.flags(ios_base::dec | ios_base::skipws);
3042 result_type __a;
3043 result_type __b;
3044 __is >> __a >> __b;
3045 if (!__is.fail())
3046 __x.param(param_type(__a, __b));
3047 return __is;
3048}
3049
Howard Hinnant30a840f2010-05-12 17:08:57 +00003050// uniform_real_distribution
3051
Howard Hinnantbc8d3f92010-05-11 19:42:16 +00003052template<class _RealType = double>
3053class uniform_real_distribution
3054{
3055public:
3056 // types
3057 typedef _RealType result_type;
3058
3059 class param_type
3060 {
3061 result_type __a_;
3062 result_type __b_;
3063 public:
3064 typedef uniform_real_distribution distribution_type;
3065
3066 explicit param_type(result_type __a = 0,
3067 result_type __b = 1)
3068 : __a_(__a), __b_(__b) {}
3069
3070 result_type a() const {return __a_;}
3071 result_type b() const {return __b_;}
3072
3073 friend bool operator==(const param_type& __x, const param_type& __y)
3074 {return __x.__a_ == __y.__a_ && __x.__b_ == __y.__b_;}
3075 friend bool operator!=(const param_type& __x, const param_type& __y)
3076 {return !(__x == __y);}
3077 };
3078
3079private:
3080 param_type __p_;
3081
3082public:
3083 // constructors and reset functions
3084 explicit uniform_real_distribution(result_type __a = 0, result_type __b = 1)
3085 : __p_(param_type(__a, __b)) {}
3086 explicit uniform_real_distribution(const param_type& __p) : __p_(__p) {}
3087 void reset() {}
3088
3089 // generating functions
3090 template<class _URNG> result_type operator()(_URNG& __g)
3091 {return (*this)(__g, __p_);}
3092 template<class _URNG> result_type operator()(_URNG& __g, const param_type& __p);
3093
3094 // property functions
3095 result_type a() const {return __p_.a();}
3096 result_type b() const {return __p_.b();}
3097
3098 param_type param() const {return __p_;}
3099 void param(const param_type& __p) {__p_ = __p;}
3100
3101 result_type min() const {return a();}
3102 result_type max() const {return b();}
3103
3104 friend bool operator==(const uniform_real_distribution& __x,
3105 const uniform_real_distribution& __y)
3106 {return __x.__p_ == __y.__p_;}
3107 friend bool operator!=(const uniform_real_distribution& __x,
3108 const uniform_real_distribution& __y)
3109 {return !(__x == __y);}
3110};
3111
3112template<class _RealType>
3113template<class _URNG>
3114inline
3115typename uniform_real_distribution<_RealType>::result_type
3116uniform_real_distribution<_RealType>::operator()(_URNG& __g, const param_type& __p)
3117{
3118 return (__p.b() - __p.a())
3119 * _STD::generate_canonical<_RealType, numeric_limits<_RealType>::digits>(__g)
3120 + __p.a();
3121}
3122
3123template <class _CharT, class _Traits, class _RT>
3124basic_ostream<_CharT, _Traits>&
3125operator<<(basic_ostream<_CharT, _Traits>& __os,
3126 const uniform_real_distribution<_RT>& __x)
3127{
3128 __save_flags<_CharT, _Traits> _(__os);
3129 __os.flags(ios_base::dec | ios_base::left);
3130 _CharT __sp = __os.widen(' ');
3131 __os.fill(__sp);
3132 return __os << __x.a() << __sp << __x.b();
3133}
3134
3135template <class _CharT, class _Traits, class _RT>
3136basic_istream<_CharT, _Traits>&
3137operator>>(basic_istream<_CharT, _Traits>& __is,
3138 uniform_real_distribution<_RT>& __x)
3139{
3140 typedef uniform_real_distribution<_RT> _Eng;
3141 typedef typename _Eng::result_type result_type;
3142 typedef typename _Eng::param_type param_type;
3143 __save_flags<_CharT, _Traits> _(__is);
3144 __is.flags(ios_base::dec | ios_base::skipws);
3145 result_type __a;
3146 result_type __b;
3147 __is >> __a >> __b;
3148 if (!__is.fail())
3149 __x.param(param_type(__a, __b));
3150 return __is;
3151}
3152
Howard Hinnant30a840f2010-05-12 17:08:57 +00003153// bernoulli_distribution
3154
Howard Hinnantbc8d3f92010-05-11 19:42:16 +00003155class bernoulli_distribution
3156{
3157public:
3158 // types
3159 typedef bool result_type;
3160
3161 class param_type
3162 {
3163 double __p_;
3164 public:
3165 typedef bernoulli_distribution distribution_type;
3166
3167 explicit param_type(double __p = 0.5) : __p_(__p) {}
3168
3169 double p() const {return __p_;}
3170
3171 friend bool operator==(const param_type& __x, const param_type& __y)
3172 {return __x.__p_ == __y.__p_;}
3173 friend bool operator!=(const param_type& __x, const param_type& __y)
3174 {return !(__x == __y);}
3175 };
3176
3177private:
3178 param_type __p_;
3179
3180public:
3181 // constructors and reset functions
3182 explicit bernoulli_distribution(double __p = 0.5)
3183 : __p_(param_type(__p)) {}
3184 explicit bernoulli_distribution(const param_type& __p) : __p_(__p) {}
3185 void reset() {}
3186
3187 // generating functions
3188 template<class _URNG> result_type operator()(_URNG& __g)
3189 {return (*this)(__g, __p_);}
Howard Hinnant03aad812010-05-11 23:26:59 +00003190 template<class _URNG> result_type operator()(_URNG& __g, const param_type& __p);
Howard Hinnantbc8d3f92010-05-11 19:42:16 +00003191
3192 // property functions
3193 double p() const {return __p_.p();}
3194
3195 param_type param() const {return __p_;}
3196 void param(const param_type& __p) {__p_ = __p;}
3197
3198 result_type min() const {return false;}
3199 result_type max() const {return true;}
3200
3201 friend bool operator==(const bernoulli_distribution& __x,
3202 const bernoulli_distribution& __y)
3203 {return __x.__p_ == __y.__p_;}
3204 friend bool operator!=(const bernoulli_distribution& __x,
3205 const bernoulli_distribution& __y)
3206 {return !(__x == __y);}
3207};
3208
Howard Hinnant03aad812010-05-11 23:26:59 +00003209template<class _URNG>
3210inline
3211bernoulli_distribution::result_type
3212bernoulli_distribution::operator()(_URNG& __g, const param_type& __p)
3213{
3214 return (__g() - __g.min()) < __p.p() * (__g.max() - __g.min() + 1.);
3215}
3216
Howard Hinnantbc8d3f92010-05-11 19:42:16 +00003217template <class _CharT, class _Traits>
3218basic_ostream<_CharT, _Traits>&
3219operator<<(basic_ostream<_CharT, _Traits>& __os, const bernoulli_distribution& __x)
3220{
3221 __save_flags<_CharT, _Traits> _(__os);
3222 __os.flags(ios_base::dec | ios_base::left);
3223 _CharT __sp = __os.widen(' ');
3224 __os.fill(__sp);
3225 return __os << __x.p();
3226}
3227
3228template <class _CharT, class _Traits>
3229basic_istream<_CharT, _Traits>&
3230operator>>(basic_istream<_CharT, _Traits>& __is, bernoulli_distribution& __x)
3231{
3232 typedef bernoulli_distribution _Eng;
3233 typedef typename _Eng::param_type param_type;
3234 __save_flags<_CharT, _Traits> _(__is);
3235 __is.flags(ios_base::dec | ios_base::skipws);
3236 double __p;
3237 __is >> __p;
3238 if (!__is.fail())
3239 __x.param(param_type(__p));
3240 return __is;
3241}
3242
Howard Hinnant30a840f2010-05-12 17:08:57 +00003243// binomial_distribution
3244
Howard Hinnant03aad812010-05-11 23:26:59 +00003245template<class _IntType = int>
3246class binomial_distribution
3247{
3248public:
3249 // types
3250 typedef _IntType result_type;
3251
3252 class param_type
3253 {
3254 result_type __t_;
3255 double __p_;
Howard Hinnant6add8dd2010-05-15 21:36:23 +00003256 double __pr_;
3257 double __odds_ratio_;
3258 result_type __r0_;
Howard Hinnant03aad812010-05-11 23:26:59 +00003259 public:
3260 typedef binomial_distribution distribution_type;
3261
Howard Hinnant6add8dd2010-05-15 21:36:23 +00003262 explicit param_type(result_type __t = 1, double __p = 0.5);
Howard Hinnant03aad812010-05-11 23:26:59 +00003263
3264 result_type t() const {return __t_;}
3265 double p() const {return __p_;}
3266
3267 friend bool operator==(const param_type& __x, const param_type& __y)
3268 {return __x.__t_ == __y.__t_ && __x.__p_ == __y.__p_;}
3269 friend bool operator!=(const param_type& __x, const param_type& __y)
3270 {return !(__x == __y);}
Howard Hinnant6add8dd2010-05-15 21:36:23 +00003271
3272 friend class binomial_distribution;
Howard Hinnant03aad812010-05-11 23:26:59 +00003273 };
3274
3275private:
3276 param_type __p_;
3277
3278public:
3279 // constructors and reset functions
3280 explicit binomial_distribution(result_type __t = 1, double __p = 0.5)
3281 : __p_(param_type(__t, __p)) {}
3282 explicit binomial_distribution(const param_type& __p) : __p_(__p) {}
3283 void reset() {}
3284
3285 // generating functions
3286 template<class _URNG> result_type operator()(_URNG& __g)
3287 {return (*this)(__g, __p_);}
3288 template<class _URNG> result_type operator()(_URNG& __g, const param_type& __p);
3289
3290 // property functions
3291 result_type t() const {return __p_.t();}
3292 double p() const {return __p_.p();}
3293
3294 param_type param() const {return __p_;}
3295 void param(const param_type& __p) {__p_ = __p;}
3296
3297 result_type min() const {return 0;}
3298 result_type max() const {return t();}
3299
3300 friend bool operator==(const binomial_distribution& __x,
3301 const binomial_distribution& __y)
3302 {return __x.__p_ == __y.__p_;}
3303 friend bool operator!=(const binomial_distribution& __x,
3304 const binomial_distribution& __y)
3305 {return !(__x == __y);}
3306};
3307
3308template<class _IntType>
Howard Hinnant6add8dd2010-05-15 21:36:23 +00003309binomial_distribution<_IntType>::param_type::param_type(result_type __t, double __p)
3310 : __t_(__t), __p_(__p)
3311{
3312 if (0 < __p_ && __p_ < 1)
3313 {
3314 __r0_ = static_cast<result_type>((__t_ + 1) * __p_);
3315 __pr_ = _STD::exp(_STD::lgamma(__t_ + 1.) - _STD::lgamma(__r0_ + 1.) -
3316 _STD::lgamma(__t_ - __r0_ + 1.) + __r0_ * _STD::log(__p_) +
3317 (__t_ - __r0_) * _STD::log(1 - __p_));
3318 __odds_ratio_ = __p_ / (1 - __p_);
3319 }
3320}
3321
3322template<class _IntType>
Howard Hinnant03aad812010-05-11 23:26:59 +00003323template<class _URNG>
3324_IntType
Howard Hinnant6add8dd2010-05-15 21:36:23 +00003325binomial_distribution<_IntType>::operator()(_URNG& __g, const param_type& __pr)
Howard Hinnant03aad812010-05-11 23:26:59 +00003326{
Howard Hinnant6add8dd2010-05-15 21:36:23 +00003327 if (__pr.__t_ == 0 || __pr.__p_ == 0)
3328 return 0;
3329 if (__pr.__p_ == 1)
3330 return __pr.__t_;
3331 uniform_real_distribution<double> __gen;
3332 double __u = __gen(__g) - __pr.__pr_;
3333 if (__u < 0)
3334 return __pr.__r0_;
3335 double __pu = __pr.__pr_;
3336 double __pd = __pu;
3337 result_type __ru = __pr.__r0_;
3338 result_type __rd = __ru;
3339 while (true)
3340 {
3341 if (__rd >= 1)
3342 {
3343 __pd *= __rd / (__pr.__odds_ratio_ * (__pr.__t_ - __rd + 1));
3344 __u -= __pd;
3345 if (__u < 0)
3346 return __rd - 1;
3347 }
3348 --__rd;
3349 ++__ru;
3350 if (__ru <= __pr.__t_)
3351 {
3352 __pu *= (__pr.__t_ - __ru + 1) * __pr.__odds_ratio_ / __ru;
3353 __u -= __pu;
3354 if (__u < 0)
3355 return __ru;
3356 }
3357 }
Howard Hinnant03aad812010-05-11 23:26:59 +00003358}
3359
3360template <class _CharT, class _Traits, class _IntType>
3361basic_ostream<_CharT, _Traits>&
3362operator<<(basic_ostream<_CharT, _Traits>& __os,
3363 const binomial_distribution<_IntType>& __x)
3364{
3365 __save_flags<_CharT, _Traits> _(__os);
3366 __os.flags(ios_base::dec | ios_base::left);
3367 _CharT __sp = __os.widen(' ');
3368 __os.fill(__sp);
3369 return __os << __x.t() << __sp << __x.p();
3370}
3371
3372template <class _CharT, class _Traits, class _IntType>
3373basic_istream<_CharT, _Traits>&
3374operator>>(basic_istream<_CharT, _Traits>& __is,
3375 binomial_distribution<_IntType>& __x)
3376{
3377 typedef binomial_distribution<_IntType> _Eng;
3378 typedef typename _Eng::result_type result_type;
3379 typedef typename _Eng::param_type param_type;
3380 __save_flags<_CharT, _Traits> _(__is);
3381 __is.flags(ios_base::dec | ios_base::skipws);
3382 result_type __t;
3383 double __p;
3384 __is >> __t >> __p;
3385 if (!__is.fail())
3386 __x.param(param_type(__t, __p));
3387 return __is;
3388}
3389
Howard Hinnant30a840f2010-05-12 17:08:57 +00003390// exponential_distribution
3391
3392template<class _RealType = double>
3393class exponential_distribution
3394{
3395public:
3396 // types
3397 typedef _RealType result_type;
3398
3399 class param_type
3400 {
3401 result_type __lambda_;
3402 public:
3403 typedef exponential_distribution distribution_type;
3404
3405 explicit param_type(result_type __lambda = 1) : __lambda_(__lambda) {}
3406
3407 result_type lambda() const {return __lambda_;}
3408
3409 friend bool operator==(const param_type& __x, const param_type& __y)
3410 {return __x.__lambda_ == __y.__lambda_;}
3411 friend bool operator!=(const param_type& __x, const param_type& __y)
3412 {return !(__x == __y);}
3413 };
3414
3415private:
3416 param_type __p_;
3417
3418public:
3419 // constructors and reset functions
3420 explicit exponential_distribution(result_type __lambda = 1)
3421 : __p_(param_type(__lambda)) {}
3422 explicit exponential_distribution(const param_type& __p) : __p_(__p) {}
3423 void reset() {}
3424
3425 // generating functions
3426 template<class _URNG> result_type operator()(_URNG& __g)
3427 {return (*this)(__g, __p_);}
3428 template<class _URNG> result_type operator()(_URNG& __g, const param_type& __p);
3429
3430 // property functions
3431 result_type lambda() const {return __p_.lambda();}
3432
3433 param_type param() const {return __p_;}
3434 void param(const param_type& __p) {__p_ = __p;}
3435
3436 result_type min() const {return 0;}
3437 result_type max() const
3438 {return -std::log(1-std::nextafter(result_type(1), result_type(-1))) /
3439 __p_.lambda();}
3440
3441 friend bool operator==(const exponential_distribution& __x,
3442 const exponential_distribution& __y)
3443 {return __x.__p_ == __y.__p_;}
3444 friend bool operator!=(const exponential_distribution& __x,
3445 const exponential_distribution& __y)
3446 {return !(__x == __y);}
3447};
3448
3449template <class _RealType>
3450template<class _URNG>
3451_RealType
3452exponential_distribution<_RealType>::operator()(_URNG& __g, const param_type& __p)
3453{
3454 return -_STD::log
3455 (
3456 result_type(1) -
3457 _STD::generate_canonical<result_type,
3458 numeric_limits<result_type>::digits>(__g)
3459 )
3460 / __p.lambda();
3461}
3462
3463template <class _CharT, class _Traits, class _RealType>
3464basic_ostream<_CharT, _Traits>&
3465operator<<(basic_ostream<_CharT, _Traits>& __os,
3466 const exponential_distribution<_RealType>& __x)
3467{
3468 __save_flags<_CharT, _Traits> _(__os);
3469 __os.flags(ios_base::dec | ios_base::left);
3470 return __os << __x.lambda();
3471}
3472
3473template <class _CharT, class _Traits, class _RealType>
3474basic_istream<_CharT, _Traits>&
3475operator>>(basic_istream<_CharT, _Traits>& __is,
3476 exponential_distribution<_RealType>& __x)
3477{
3478 typedef exponential_distribution<_RealType> _Eng;
3479 typedef typename _Eng::result_type result_type;
3480 typedef typename _Eng::param_type param_type;
3481 __save_flags<_CharT, _Traits> _(__is);
3482 __is.flags(ios_base::dec | ios_base::skipws);
3483 result_type __lambda;
3484 __is >> __lambda;
3485 if (!__is.fail())
3486 __x.param(param_type(__lambda));
3487 return __is;
3488}
3489
Howard Hinnant6add8dd2010-05-15 21:36:23 +00003490// normal_distribution
3491
3492template<class _RealType = double>
3493class normal_distribution
3494{
3495public:
3496 // types
3497 typedef _RealType result_type;
3498
3499 class param_type
3500 {
3501 result_type __mean_;
3502 result_type __stddev_;
3503 public:
3504 typedef normal_distribution distribution_type;
3505
3506 explicit param_type(result_type __mean = 0, result_type __stddev = 1)
3507 : __mean_(__mean), __stddev_(__stddev) {}
3508
3509 result_type mean() const {return __mean_;}
3510 result_type stddev() const {return __stddev_;}
3511
3512 friend bool operator==(const param_type& __x, const param_type& __y)
3513 {return __x.__mean_ == __y.__mean_ && __x.__stddev_ == __y.__stddev_;}
3514 friend bool operator!=(const param_type& __x, const param_type& __y)
3515 {return !(__x == __y);}
3516 };
3517
3518private:
3519 param_type __p_;
3520 result_type _V_;
3521 bool _V_hot_;
3522
3523public:
3524 // constructors and reset functions
3525 explicit normal_distribution(result_type __mean = 0, result_type __stddev = 1)
3526 : __p_(param_type(__mean, __stddev)), _V_hot_(false) {}
3527 explicit normal_distribution(const param_type& __p)
3528 : __p_(__p), _V_hot_(false) {}
3529 void reset() {_V_hot_ = false;}
3530
3531 // generating functions
3532 template<class _URNG> result_type operator()(_URNG& __g)
3533 {return (*this)(__g, __p_);}
3534 template<class _URNG> result_type operator()(_URNG& __g, const param_type& __p);
3535
3536 // property functions
3537 result_type mean() const {return __p_.mean();}
3538 result_type stddev() const {return __p_.stddev();}
3539
3540 param_type param() const {return __p_;}
3541 void param(const param_type& __p) {__p_ = __p;}
3542
3543 result_type min() const {return -numeric_limits<result_type>::infinity();}
3544 result_type max() const {return numeric_limits<result_type>::infinity();}
3545
3546 friend bool operator==(const normal_distribution& __x,
3547 const normal_distribution& __y)
3548 {return __x.__p_ == __y.__p_ && __x._V_hot_ == __y._V_hot_ &&
3549 (!__x._V_hot_ || __x._V_ == __y._V_);}
3550 friend bool operator!=(const normal_distribution& __x,
3551 const normal_distribution& __y)
3552 {return !(__x == __y);}
3553
3554 template <class _CharT, class _Traits, class _RT>
3555 friend
3556 basic_ostream<_CharT, _Traits>&
3557 operator<<(basic_ostream<_CharT, _Traits>& __os,
3558 const normal_distribution<_RT>& __x);
3559
3560 template <class _CharT, class _Traits, class _RT>
3561 friend
3562 basic_istream<_CharT, _Traits>&
3563 operator>>(basic_istream<_CharT, _Traits>& __is,
3564 normal_distribution<_RT>& __x);
3565};
3566
3567template <class _RealType>
3568template<class _URNG>
3569_RealType
3570normal_distribution<_RealType>::operator()(_URNG& __g, const param_type& __p)
3571{
3572 result_type _U;
3573 if (_V_hot_)
3574 {
3575 _V_hot_ = false;
3576 _U = _V_;
3577 }
3578 else
3579 {
3580 uniform_real_distribution<result_type> _Uni(-1, 1);
3581 result_type __u;
3582 result_type __v;
3583 result_type __s;
3584 do
3585 {
3586 __u = _Uni(__g);
3587 __v = _Uni(__g);
3588 __s = __u * __u + __v * __v;
3589 } while (__s > 1 || __s == 0);
3590 result_type _F = _STD::sqrt(-2 * _STD::log(__s) / __s);
3591 _V_ = __v * _F;
3592 _V_hot_ = true;
3593 _U = __u * _F;
3594 }
3595 return _U * __p.stddev() + __p.mean();
3596}
3597
3598template <class _CharT, class _Traits, class _RT>
3599basic_ostream<_CharT, _Traits>&
3600operator<<(basic_ostream<_CharT, _Traits>& __os,
3601 const normal_distribution<_RT>& __x)
3602{
3603 __save_flags<_CharT, _Traits> _(__os);
3604 __os.flags(ios_base::dec | ios_base::left);
3605 _CharT __sp = __os.widen(' ');
3606 __os.fill(__sp);
3607 __os << __x.mean() << __sp << __x.stddev() << __sp << __x._V_hot_;
3608 if (__x._V_hot_)
3609 __os << __sp << __x._V_;
3610 return __os;
3611}
3612
3613template <class _CharT, class _Traits, class _RT>
3614basic_istream<_CharT, _Traits>&
3615operator>>(basic_istream<_CharT, _Traits>& __is,
3616 normal_distribution<_RT>& __x)
3617{
3618 typedef normal_distribution<_RT> _Eng;
3619 typedef typename _Eng::result_type result_type;
3620 typedef typename _Eng::param_type param_type;
3621 __save_flags<_CharT, _Traits> _(__is);
3622 __is.flags(ios_base::dec | ios_base::skipws);
3623 result_type __mean;
3624 result_type __stddev;
3625 result_type _V = 0;
3626 bool _V_hot = false;
3627 __is >> __mean >> __stddev >> _V_hot;
3628 if (_V_hot)
3629 __is >> _V;
3630 if (!__is.fail())
3631 {
3632 __x.param(param_type(__mean, __stddev));
3633 __x._V_hot_ = _V_hot;
3634 __x._V_ = _V;
3635 }
3636 return __is;
3637}
3638
3639// poisson_distribution
3640
3641template<class _IntType = int>
3642class poisson_distribution
3643{
3644public:
3645 // types
3646 typedef _IntType result_type;
3647
3648 class param_type
3649 {
3650 double __mean_;
3651 double __s_;
3652 double __d_;
3653 double __l_;
3654 double __omega_;
3655 double __c0_;
3656 double __c1_;
3657 double __c2_;
3658 double __c3_;
3659 double __c_;
3660
3661 public:
3662 typedef poisson_distribution distribution_type;
3663
3664 explicit param_type(double __mean = 1.0);
3665
3666 double mean() const {return __mean_;}
3667
3668 friend bool operator==(const param_type& __x, const param_type& __y)
3669 {return __x.__mean_ == __y.__mean_;}
3670 friend bool operator!=(const param_type& __x, const param_type& __y)
3671 {return !(__x == __y);}
3672
3673 friend class poisson_distribution;
3674 };
3675
3676private:
3677 param_type __p_;
3678
3679public:
3680 // constructors and reset functions
3681 explicit poisson_distribution(double __mean = 1.0) : __p_(__mean) {}
3682 explicit poisson_distribution(const param_type& __p) : __p_(__p) {}
3683 void reset() {}
3684
3685 // generating functions
3686 template<class _URNG> result_type operator()(_URNG& __g)
3687 {return (*this)(__g, __p_);}
3688 template<class _URNG> result_type operator()(_URNG& __g, const param_type& __p);
3689
3690 // property functions
3691 double mean() const {return __p_.mean();}
3692
3693 param_type param() const {return __p_;}
3694 void param(const param_type& __p) {__p_ = __p;}
3695
3696 result_type min() const {return 0;}
3697 result_type max() const {return numeric_limits<result_type>::max();}
3698
3699 friend bool operator==(const poisson_distribution& __x,
3700 const poisson_distribution& __y)
3701 {return __x.__p_ == __y.__p_;}
3702 friend bool operator!=(const poisson_distribution& __x,
3703 const poisson_distribution& __y)
3704 {return !(__x == __y);}
3705};
3706
3707template<class _IntType>
3708poisson_distribution<_IntType>::param_type::param_type(double __mean)
3709 : __mean_(__mean)
3710{
3711 if (__mean_ < 10)
3712 {
3713 __s_ = 0;
3714 __d_ = 0;
3715 __l_ = _STD::exp(-__mean_);
3716 __omega_ = 0;
3717 __c3_ = 0;
3718 __c2_ = 0;
3719 __c1_ = 0;
3720 __c0_ = 0;
3721 __c_ = 0;
3722 }
3723 else
3724 {
3725 __s_ = _STD::sqrt(__mean_);
3726 __d_ = 6 * __mean_ * __mean_;
3727 __l_ = static_cast<result_type>(__mean_ - 1.1484);
3728 __omega_ = .3989423 / __s_;
3729 double __b1_ = .4166667E-1 / __mean_;
3730 double __b2_ = .3 * __b1_ * __b1_;
3731 __c3_ = .1428571 * __b1_ * __b2_;
3732 __c2_ = __b2_ - 15. * __c3_;
3733 __c1_ = __b1_ - 6. * __b2_ + 45. * __c3_;
3734 __c0_ = 1. - __b1_ + 3. * __b2_ - 15. * __c3_;
3735 __c_ = .1069 / __mean_;
3736 }
3737}
3738
3739template <class _IntType>
3740template<class _URNG>
3741_IntType
3742poisson_distribution<_IntType>::operator()(_URNG& __urng, const param_type& __pr)
3743{
3744 result_type __x;
3745 uniform_real_distribution<double> __urd;
3746 if (__pr.__mean_ <= 10)
3747 {
3748 __x = 0;
3749 for (double __p = __urd(__urng); __p > __pr.__l_; ++__x)
3750 __p *= __urd(__urng);
3751 }
3752 else
3753 {
3754 double __difmuk;
3755 double __g = __pr.__mean_ + __pr.__s_ * normal_distribution<double>()(__urng);
3756 double __u;
3757 if (__g > 0)
3758 {
3759 __x = static_cast<result_type>(__g);
3760 if (__x >= __pr.__l_)
3761 return __x;
3762 __difmuk = __pr.__mean_ - __x;
3763 __u = __urd(__urng);
3764 if (__pr.__d_ * __u >= __difmuk * __difmuk * __difmuk)
3765 return __x;
3766 }
3767 exponential_distribution<double> __edist;
3768 for (bool __using_exp_dist = false; true; __using_exp_dist = true)
3769 {
3770 double __e;
3771 if (__using_exp_dist || __g < 0)
3772 {
3773 double __t;
3774 do
3775 {
3776 __e = __edist(__urng);
3777 __u = __urd(__urng);
3778 __u += __u - 1;
3779 __t = 1.8 + (__u < 0 ? -__e : __e);
3780 } while (__t <= -.6744);
3781 __x = __pr.__mean_ + __pr.__s_ * __t;
3782 __difmuk = __pr.__mean_ - __x;
3783 __using_exp_dist = true;
3784 }
3785 double __px;
3786 double __py;
3787 if (__x < 10)
3788 {
3789 const result_type __fac[] = {1, 1, 2, 6, 24, 120, 720, 5040,
3790 40320, 362880};
3791 __px = -__pr.__mean_;
3792 __py = _STD::pow(__pr.__mean_, (double)__x) / __fac[__x];
3793 }
3794 else
3795 {
3796 double __del = .8333333E-1 / __x;
3797 __del -= 4.8 * __del * __del * __del;
3798 double __v = __difmuk / __x;
3799 if (_STD::abs(__v) > 0.25)
3800 __px = __x * _STD::log(1 + __v) - __difmuk - __del;
3801 else
3802 __px = __x * __v * __v * (((((((.1250060 * __v + -.1384794) *
3803 __v + .1421878) * __v + -.1661269) * __v + .2000118) *
3804 __v + -.2500068) * __v + .3333333) * __v + -.5) - __del;
3805 __py = .3989423 / _STD::sqrt(__x);
3806 }
3807 double __r = (0.5 - __difmuk) / __pr.__s_;
3808 double __r2 = __r * __r;
3809 double __fx = -0.5 * __r2;
3810 double __fy = __pr.__omega_ * (((__pr.__c3_ * __r2 + __pr.__c2_) *
3811 __r2 + __pr.__c1_) * __r2 + __pr.__c0_);
3812 if (__using_exp_dist)
3813 {
3814 if (__pr.__c_ * _STD::abs(__u) <= __py * _STD::exp(__px + __e) -
3815 __fy * _STD::exp(__fx + __e))
3816 break;
3817 }
3818 else
3819 {
3820 if (__fy - __u * __fy <= __py * _STD::exp(__px - __fx))
3821 break;
3822 }
3823 }
3824 }
3825 return __x;
3826}
3827
3828template <class _CharT, class _Traits, class _IntType>
3829basic_ostream<_CharT, _Traits>&
3830operator<<(basic_ostream<_CharT, _Traits>& __os,
3831 const poisson_distribution<_IntType>& __x)
3832{
3833 __save_flags<_CharT, _Traits> _(__os);
3834 __os.flags(ios_base::dec | ios_base::left);
3835 return __os << __x.mean();
3836}
3837
3838template <class _CharT, class _Traits, class _IntType>
3839basic_istream<_CharT, _Traits>&
3840operator>>(basic_istream<_CharT, _Traits>& __is,
3841 poisson_distribution<_IntType>& __x)
3842{
3843 typedef poisson_distribution<_IntType> _Eng;
3844 typedef typename _Eng::param_type param_type;
3845 __save_flags<_CharT, _Traits> _(__is);
3846 __is.flags(ios_base::dec | ios_base::skipws);
3847 double __mean;
3848 __is >> __mean;
3849 if (!__is.fail())
3850 __x.param(param_type(__mean));
3851 return __is;
3852}
3853
Howard Hinnant9de6e302010-05-16 01:09:02 +00003854// weibull_distribution
3855
3856template<class _RealType = double>
3857class weibull_distribution
3858{
3859public:
3860 // types
3861 typedef _RealType result_type;
3862
3863 class param_type
3864 {
3865 result_type __a_;
3866 result_type __b_;
3867 public:
3868 typedef weibull_distribution distribution_type;
3869
3870 explicit param_type(result_type __a = 1, result_type __b = 1)
3871 : __a_(__a), __b_(__b) {}
3872
3873 result_type a() const {return __a_;}
3874 result_type b() const {return __b_;}
3875
3876 friend bool operator==(const param_type& __x, const param_type& __y)
3877 {return __x.__a_ == __y.__a_ && __x.__b_ == __y.__b_;}
3878 friend bool operator!=(const param_type& __x, const param_type& __y)
3879 {return !(__x == __y);}
3880 };
3881
3882private:
3883 param_type __p_;
3884
3885public:
3886 // constructor and reset functions
3887 explicit weibull_distribution(result_type __a = 1, result_type __b = 1)
3888 : __p_(param_type(__a, __b)) {}
3889 explicit weibull_distribution(const param_type& __p)
3890 : __p_(__p) {}
3891 void reset() {}
3892
3893 // generating functions
3894 template<class _URNG> result_type operator()(_URNG& __g)
3895 {return (*this)(__g, __p_);}
3896 template<class _URNG> result_type operator()(_URNG& __g, const param_type& __p)
3897 {return __p.b() *
3898 _STD::pow(exponential_distribution<result_type>()(__g), 1/__p.a());}
3899
3900 // property functions
3901 result_type a() const {return __p_.a();}
3902 result_type b() const {return __p_.b();}
3903
3904 param_type param() const {return __p_;}
3905 void param(const param_type& __p) {__p_ = __p;}
3906
3907 result_type min() const {return 0;}
3908 result_type max() const {return numeric_limits<result_type>::infinity();}
3909
3910
3911 friend bool operator==(const weibull_distribution& __x,
3912 const weibull_distribution& __y)
3913 {return __x.__p_ == __y.__p_;}
3914 friend bool operator!=(const weibull_distribution& __x,
3915 const weibull_distribution& __y)
3916 {return !(__x == __y);}
3917};
3918
3919template <class _CharT, class _Traits, class _RT>
3920basic_ostream<_CharT, _Traits>&
3921operator<<(basic_ostream<_CharT, _Traits>& __os,
3922 const weibull_distribution<_RT>& __x)
3923{
3924 __save_flags<_CharT, _Traits> _(__os);
3925 __os.flags(ios_base::dec | ios_base::left);
3926 _CharT __sp = __os.widen(' ');
3927 __os.fill(__sp);
3928 __os << __x.a() << __sp << __x.b();
3929 return __os;
3930}
3931
3932template <class _CharT, class _Traits, class _RT>
3933basic_istream<_CharT, _Traits>&
3934operator>>(basic_istream<_CharT, _Traits>& __is,
3935 weibull_distribution<_RT>& __x)
3936{
3937 typedef weibull_distribution<_RT> _Eng;
3938 typedef typename _Eng::result_type result_type;
3939 typedef typename _Eng::param_type param_type;
3940 __save_flags<_CharT, _Traits> _(__is);
3941 __is.flags(ios_base::dec | ios_base::skipws);
3942 result_type __a;
3943 result_type __b;
3944 __is >> __a >> __b;
3945 if (!__is.fail())
3946 __x.param(param_type(__a, __b));
3947 return __is;
3948}
3949
Howard Hinnantc7c49132010-05-13 17:58:28 +00003950// gamma_distribution
3951
3952template<class _RealType = double>
3953class gamma_distribution
3954{
3955public:
3956 // types
3957 typedef _RealType result_type;
3958
3959 class param_type
3960 {
3961 result_type __alpha_;
3962 result_type __beta_;
3963 public:
3964 typedef gamma_distribution distribution_type;
3965
3966 explicit param_type(result_type __alpha = 1, result_type __beta = 1)
3967 : __alpha_(__alpha), __beta_(__beta) {}
3968
3969 result_type alpha() const {return __alpha_;}
3970 result_type beta() const {return __beta_;}
3971
3972 friend bool operator==(const param_type& __x, const param_type& __y)
3973 {return __x.__alpha_ == __y.__alpha_ && __x.__beta_ == __y.__beta_;}
3974 friend bool operator!=(const param_type& __x, const param_type& __y)
3975 {return !(__x == __y);}
3976 };
3977
3978private:
3979 param_type __p_;
3980
3981public:
3982 // constructors and reset functions
3983 explicit gamma_distribution(result_type __alpha = 1, result_type __beta = 1)
3984 : __p_(param_type(__alpha, __beta)) {}
3985 explicit gamma_distribution(const param_type& __p)
3986 : __p_(__p) {}
3987 void reset() {}
3988
3989 // generating functions
3990 template<class _URNG> result_type operator()(_URNG& __g)
3991 {return (*this)(__g, __p_);}
3992 template<class _URNG> result_type operator()(_URNG& __g, const param_type& __p);
3993
3994 // property functions
3995 result_type alpha() const {return __p_.alpha();}
3996 result_type beta() const {return __p_.beta();}
3997
3998 param_type param() const {return __p_;}
3999 void param(const param_type& __p) {__p_ = __p;}
4000
4001 result_type min() const {return 0;}
4002 result_type max() const {return numeric_limits<result_type>::infinity();}
4003
4004 friend bool operator==(const gamma_distribution& __x,
4005 const gamma_distribution& __y)
4006 {return __x.__p_ == __y.__p_;}
4007 friend bool operator!=(const gamma_distribution& __x,
4008 const gamma_distribution& __y)
4009 {return !(__x == __y);}
4010};
4011
4012template <class _RealType>
4013template<class _URNG>
4014_RealType
4015gamma_distribution<_RealType>::operator()(_URNG& __g, const param_type& __p)
4016{
Howard Hinnantf417abe2010-05-14 18:43:10 +00004017 result_type __a = __p.alpha();
4018 uniform_real_distribution<result_type> __gen(0, 1);
4019 exponential_distribution<result_type> __egen;
4020 result_type __x;
Howard Hinnantc7c49132010-05-13 17:58:28 +00004021 if (__a == 1)
Howard Hinnantf417abe2010-05-14 18:43:10 +00004022 __x = __egen(__g);
Howard Hinnantc7c49132010-05-13 17:58:28 +00004023 else if (__a > 1)
4024 {
4025 const result_type __b = __a - 1;
4026 const result_type __c = 3 * __a - result_type(0.75);
Howard Hinnantc7c49132010-05-13 17:58:28 +00004027 while (true)
4028 {
4029 const result_type __u = __gen(__g);
4030 const result_type __v = __gen(__g);
4031 const result_type __w = __u * (1 - __u);
Howard Hinnantf417abe2010-05-14 18:43:10 +00004032 if (__w != 0)
Howard Hinnantc7c49132010-05-13 17:58:28 +00004033 {
4034 const result_type __y = _STD::sqrt(__c / __w) *
4035 (__u - result_type(0.5));
4036 __x = __b + __y;
4037 if (__x >= 0)
4038 {
4039 const result_type __z = 64 * __w * __w * __w * __v * __v;
4040 if (__z <= 1 - 2 * __y * __y / __x)
4041 break;
4042 if (_STD::log(__z) <= 2 * (__b * _STD::log(__x / __b) - __y))
4043 break;
4044 }
4045 }
4046 }
Howard Hinnantc7c49132010-05-13 17:58:28 +00004047 }
Howard Hinnantf417abe2010-05-14 18:43:10 +00004048 else // __a < 1
4049 {
4050 while (true)
4051 {
4052 const result_type __u = __gen(__g);
4053 const result_type __es = __egen(__g);
4054 if (__u <= 1 - __a)
4055 {
4056 __x = _STD::pow(__u, 1 / __a);
4057 if (__x <= __es)
4058 break;
4059 }
4060 else
4061 {
4062 const result_type __e = -_STD::log((1-__u)/__a);
4063 __x = _STD::pow(1 - __a + __a * __e, 1 / __a);
4064 if (__x <= __e + __es)
4065 break;
4066 }
4067 }
4068 }
4069 return __x * __p.beta();
Howard Hinnantc7c49132010-05-13 17:58:28 +00004070}
4071
4072template <class _CharT, class _Traits, class _RT>
4073basic_ostream<_CharT, _Traits>&
4074operator<<(basic_ostream<_CharT, _Traits>& __os,
4075 const gamma_distribution<_RT>& __x)
4076{
4077 __save_flags<_CharT, _Traits> _(__os);
4078 __os.flags(ios_base::dec | ios_base::left);
4079 _CharT __sp = __os.widen(' ');
4080 __os.fill(__sp);
4081 __os << __x.alpha() << __sp << __x.beta();
4082 return __os;
4083}
4084
4085template <class _CharT, class _Traits, class _RT>
4086basic_istream<_CharT, _Traits>&
4087operator>>(basic_istream<_CharT, _Traits>& __is,
4088 gamma_distribution<_RT>& __x)
4089{
4090 typedef gamma_distribution<_RT> _Eng;
4091 typedef typename _Eng::result_type result_type;
4092 typedef typename _Eng::param_type param_type;
4093 __save_flags<_CharT, _Traits> _(__is);
4094 __is.flags(ios_base::dec | ios_base::skipws);
4095 result_type __alpha;
4096 result_type __beta;
4097 __is >> __alpha >> __beta;
4098 if (!__is.fail())
4099 __x.param(param_type(__alpha, __beta));
4100 return __is;
4101}
Howard Hinnanta64111c2010-05-12 21:02:31 +00004102
Howard Hinnant97dc2f32010-05-15 23:36:00 +00004103// chi_squared_distribution
4104
4105template<class _RealType = double>
4106class chi_squared_distribution
4107{
4108public:
4109 // types
4110 typedef _RealType result_type;
4111
4112 class param_type
4113 {
4114 result_type __n_;
4115 public:
4116 typedef chi_squared_distribution distribution_type;
4117
4118 explicit param_type(result_type __n = 1) : __n_(__n) {}
4119
4120 result_type n() const {return __n_;}
4121
4122 friend bool operator==(const param_type& __x, const param_type& __y)
4123 {return __x.__n_ == __y.__n_;}
4124 friend bool operator!=(const param_type& __x, const param_type& __y)
4125 {return !(__x == __y);}
4126 };
4127
4128private:
4129 param_type __p_;
4130
4131public:
4132 // constructor and reset functions
4133 explicit chi_squared_distribution(result_type __n = 1)
4134 : __p_(param_type(__n)) {}
4135 explicit chi_squared_distribution(const param_type& __p)
4136 : __p_(__p) {}
4137 void reset() {}
4138
4139 // generating functions
4140 template<class _URNG> result_type operator()(_URNG& __g)
4141 {return (*this)(__g, __p_);}
4142 template<class _URNG> result_type operator()(_URNG& __g, const param_type& __p)
4143 {return gamma_distribution<result_type>(__p.n() / 2, 2)(__g);}
4144
4145 // property functions
4146 result_type n() const {return __p_.n();}
4147
4148 param_type param() const {return __p_;}
4149 void param(const param_type& __p) {__p_ = __p;}
4150
4151 result_type min() const {return 0;}
4152 result_type max() const {return numeric_limits<result_type>::infinity();}
4153
4154
4155 friend bool operator==(const chi_squared_distribution& __x,
4156 const chi_squared_distribution& __y)
4157 {return __x.__p_ == __y.__p_;}
4158 friend bool operator!=(const chi_squared_distribution& __x,
4159 const chi_squared_distribution& __y)
4160 {return !(__x == __y);}
4161};
4162
4163template <class _CharT, class _Traits, class _RT>
4164basic_ostream<_CharT, _Traits>&
4165operator<<(basic_ostream<_CharT, _Traits>& __os,
4166 const chi_squared_distribution<_RT>& __x)
4167{
4168 __save_flags<_CharT, _Traits> _(__os);
4169 __os.flags(ios_base::dec | ios_base::left);
4170 __os << __x.n();
4171 return __os;
4172}
4173
4174template <class _CharT, class _Traits, class _RT>
4175basic_istream<_CharT, _Traits>&
4176operator>>(basic_istream<_CharT, _Traits>& __is,
4177 chi_squared_distribution<_RT>& __x)
4178{
4179 typedef chi_squared_distribution<_RT> _Eng;
4180 typedef typename _Eng::result_type result_type;
4181 typedef typename _Eng::param_type param_type;
4182 __save_flags<_CharT, _Traits> _(__is);
4183 __is.flags(ios_base::dec | ios_base::skipws);
4184 result_type __n;
4185 __is >> __n;
4186 if (!__is.fail())
4187 __x.param(param_type(__n));
4188 return __is;
4189}
4190
Howard Hinnantbc8d3f92010-05-11 19:42:16 +00004191_LIBCPP_END_NAMESPACE_STD
4192
4193#endif // _LIBCPP_RANDOM