simple_normal

#include 
#include 
#include 
int main() {
  unsigned s = std::random_device()(); // Set seed from random_device
  std::mt19937 g(s);                   // Initialize URNG
  std::cout << "Seed set to " << s << "\n\n";
  exrandom::unit_normal_distribution<double> N;
  std::cout
    << "Sampling doubles exactly from the unit normal distribution:\n";
  for (int i = 0; i < 48; ++i)
    std::cout << N(g) << ((i + 1) % 8 ? ' ' : '\n');
  std::cout << "\n";
}
1
2
3
4
5
6
7
8
9
10
11
12
13
14

simple_uniform

从[0,1]范围 按照均匀分布采样

#include 
#include 
#include 
int main() {
  unsigned s = std::random_device()(); // Set seed from random_device
  std::mt19937 g(s);                   // Initialize URNG
  std::cout << "Seed set to " << s << "\n\n";
  exrandom::unit_uniform_distribution<long double> N;
  std::cout
    << "Sampling long doubles exactly from the unit uniform distribution:\n";
  for (int i = 0; i < 48; ++i)
    std::cout << N(g) << ((i + 1) % 8 ? ' ' : '\n');
  std::cout << "\n";
}
1
2
3
4
5
6
7
8
9
10
11
12
13
14

simple_exponential

#include 
#include 
#include 
int main() {
  unsigned s = std::random_device()(); // Set seed from random_device
  std::mt19937 g(s);                   // Initialize URNG
  std::cout << "Seed set to " << s << "\n\n";
  exrandom::unit_exponential_distribution<float> N;
  std::cout
    << "Sampling floats exactly from the unit exponential distribution:\n";
  for (int i = 0; i < 48; ++i)
    std::cout << N(g) << ((i + 1) % 8 ? ' ' : '\n');
  std::cout << "\n";
}
1
2
3
4
5
6
7
8
9
10
11
12
13
14

simple_discrete_normal

#include 
#include 
#include 
int main() {
  unsigned s = std::random_device()(); // Set seed from random_device
  std::mt19937 g(s);                   // Initialize URNG
  std::cout << "Seed set to " << s << "\n\n";
  exrandom::discrete_normal_distribution N(1,3,129,2);
  std::cout
    << "Sampling exactly from the discrete normal distribution,\n"
    << "  mu = " << N.mu_num() << "/" << N.mu_den()
    << ", sigma = " << N.sigma_num() << "/" << N.sigma_den() << ":\n";
  for (int i = 0; i < 100; ++i)
    std::cout << N(g) << ((i + 1) % 20 ? ' ' : '\n');
  std::cout << "\n";
}
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16

sample_uniform

#include 
#include 
#include 
#include 
#include 
int main() {
  std::cout << "Sampling from the uniform distribution\n";
  unsigned s = std::random_device()(); // Set seed from random_device
  std::mt19937 g(s);                   // Initialize URNG
  std::cout << "Seed set to " << s << "\n\n";
  {
    const unsigned b = 16;
    std::cout << "Deviates rounded to floating point with base = " << b << "\n"
              << "columns are: initial u-rand, rounded float, final u-rand\n"
              << std::setprecision(std::numeric_limits<float>::digits10);
    exrandom::rand_digit<b> D;
    exrandom::unit_uniform_dist<exrandom::rand_digit<b>> N(D);
    exrandom::u_rand<exrandom::rand_digit<b>> x(N.digit_generator());
    for (int i = 0; i < 10; ++i) {
      N.generate(g, x);
      std::cout << x << " = ";
      std::cout << x.value<float>(g) << " = ";
      std::cout << x << "\n";
    }
    std::cout << "Total number of base " << b << " digits used = "
              << D.count() << "\n\n";
  }
  {
    const unsigned b = 10;
    int p = 6;
    std::cout << "Deviates rounded to fixed point with base = " << b
              << " and precision = " << p << "\n"
              << "columns are: initial u-rand, rounded fixed, final u-rand\n";
    exrandom::rand_digit<b> D;
    exrandom::unit_uniform_dist<exrandom::rand_digit<b>> N(D);
    exrandom::u_rand<exrandom::rand_digit<b>> x(N.digit_generator());
    for (int i = 0; i < 10; ++i) {
      N.generate(g, x);
      std::cout << x << " = ";
      std::cout << x.print_fixed(g,p) << " = ";
      std::cout << x << "\n";
    }
    std::cout << "Total number of base " << b << " digits used = "
              << D.count() << "\n\n";
  }
}
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46

sample_normal

#include 
#include 
#include 
#include 
#include 
int main() {
  std::cout << "Sampling from the normal distribution\n";
  unsigned s = std::random_device()(); // Set seed from random_device
  std::mt19937 g(s);                   // Initialize URNG
  std::cout << "Seed set to " << s << "\n\n";
  {
    const unsigned b = 16;
    std::cout << "Deviates rounded to floating point with base = " << b << "\n"
              << "columns are: initial u-rand, rounded float, final u-rand\n"
              << std::setprecision(std::numeric_limits<float>::digits10);
    exrandom::rand_digit<b> D;
    exrandom::unit_normal_dist<exrandom::rand_digit<b>> N(D);
    exrandom::u_rand<exrandom::rand_digit<b>> x(N.digit_generator());
    for (int i = 0; i < 10; ++i) {
      N.generate(g, x);
      std::cout << x << " = ";
      std::cout << x.value<float>(g) << " = ";
      std::cout << x << "\n";
    }
    std::cout << "Total number of base " << b << " digits used = "
              << D.count() << "\n\n";
  }
  {
    const unsigned b = 10;
    int p = 6;
    std::cout << "Deviates rounded to fixed point with base = " << b
              << " and precision = " << p << "\n"
              << "columns are: initial u-rand, rounded fixed, final u-rand\n";
    exrandom::rand_digit<b> D;
    exrandom::unit_normal_dist<exrandom::rand_digit<b>> N(D);
    exrandom::u_rand<exrandom::rand_digit<b>> x(N.digit_generator());
    for (int i = 0; i < 10; ++i) {
      N.generate(g, x);
      std::cout << x << " = ";
      std::cout << x.print_fixed(g,p) << " = ";
      std::cout << x << "\n";
    }
    std::cout << "Total number of base " << b << " digits used = "
              << D.count() << "\n\n";
  }
}
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46

sample_exponential

#include 
#include 
#include 
#include 
#include 
int main() {
  std::cout << "Sampling from the exponential distribution\n";
  unsigned s = std::random_device()(); // Set seed from random_device
  std::mt19937 g(s);                   // Initialize URNG
  std::cout << "Seed set to " << s << "\n\n";
  {
    const unsigned b = 16;
    std::cout << "Deviates rounded to floating point with base = " << b << "\n"
              << "columns are: initial u-rand, rounded float, final u-rand\n"
              << std::setprecision(std::numeric_limits<float>::digits10);
    exrandom::rand_digit<b> D;
    exrandom::unit_exponential_dist<exrandom::rand_digit<b>> N(D);
    exrandom::u_rand<exrandom::rand_digit<b>> x(N.digit_generator());
    for (int i = 0; i < 10; ++i) {
      N.generate(g, x);
      std::cout << x << " = ";
      std::cout << x.value<float>(g) << " = ";
      std::cout << x << "\n";
    }
    std::cout << "Total number of base " << b << " digits used = "
              << D.count() << "\n\n";
  }
  {
    const unsigned b = 10;
    int p = 6;
    std::cout << "Deviates rounded to fixed point with base = " << b
              << " and precision = " << p << "\n"
              << "columns are: initial u-rand, rounded fixed, final u-rand\n";
    exrandom::rand_digit<b> D;
    exrandom::unit_exponential_dist<exrandom::rand_digit<b>> N(D);
    exrandom::u_rand<exrandom::rand_digit<b>> x(N.digit_generator());
    for (int i = 0; i < 10; ++i) {
      N.generate(g, x);
      std::cout << x << " = ";
      std::cout << x.print_fixed(g,p) << " = ";
      std::cout << x << "\n";
    }
    std::cout << "Total number of base " << b << " digits used = "
              << D.count() << "\n\n";
  }
}
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46

sample_distributions

#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
template<typename Dist, typename Generator>
void histogram(long long num, Dist& d, Generator& g,
               double step, double offset) {
  std::map<int,unsigned> hist;
  const unsigned maxlength = 40; // Widest bar in histogram
  for (long long i = 0; i < num; ++i)
    ++hist[int(std::floor((double(d(g)) - offset) / step + 0.5))];
  unsigned m = 0;
  for (auto p : hist) m = (std::max)(m, p.second);
  m = (std::max)(m / maxlength, 1U);
  for (auto p : hist) {
    unsigned k = (p.second + m/2) / m;
    if (k)
      std::cout << typename Dist::result_type(p.first * step + offset) << '\t'
                << std::string(k, '*') << '\n';
  }
  std::cout << std::endl;
}
int main() {
  unsigned s = std::random_device()(); // Set seed from random_device
  std::mt19937 g(s);                   // Initialize URNG
  std::cout << "Seed set to " << s << "\n\n";
  long long num = 1000000LL;
  {
    exrandom::unit_normal_distribution<double> d;
    std::cout << "Sampling the unit normal distribution\n";
    histogram(num, d, g, 0.2, 0.0);
  }
  {
    exrandom::unit_exponential_distribution<double> d;
    std::cout << "Sampling the unit exponential distribution\n";
    histogram(num, d, g, 0.2, 0.1);
  }
  {
    exrandom::unit_uniform_distribution<double> d;
    std::cout << "Sampling the unit uniform distribution\n";
    histogram(num, d, g, 0.1, 0.05);
  }
  {
    exrandom::discrete_normal_distribution d(6, 55, 10);
    std::cout << "Sampling the discrete normal distribution"
              << " mu = " << d.mu_num() << "/" << d.mu_den()
              << " sigma = " << d.sigma_num() << "/" << d.sigma_den() << "\n";
    histogram(num, d, g, 1, 0);
  }
}
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56

sample_discrete_normal

#include 
#include 
#include 
#include 
#include 
int main() {
  std::cout << "Sampling from the discrete normal distribution\n";
  unsigned s = std::random_device()(); // Set seed from random_device
  std::mt19937 g(s);                   // Initialize URNG
  std::cout << "Seed set to " << s << "\n\n";
  const unsigned b = 2;
  exrandom::rand_digit<b> D;
  exrandom::discrete_normal_dist<exrandom::rand_digit<b>>
    N(D, 87331, 3, 16000051, 100);
  exrandom::i_rand<exrandom::rand_digit<b>> h(N.digit_generator());
  std::cout << "Parameters: "
            << " mu = " << N.mu_num() << "/" << N.mu_den()
            << " sigma = " << N.sigma_num() << "/" << N.sigma_den() << "\n"
            << "columns are: uniform range, final int\n";
  int num = 20;
  for (int i = 0; i < num; ++i) {
    N.generate(g, h);
    std::cout << h << " = ";
    std::cout << h(g) << "\n";
  }
    std::cout << "Total number of base " << b << " digits used = "
              << D.count() << "\n\n";
}
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28

tabulate_normals

#include 
#include 
#include 
#include 
#include 
#include 
#include 
int main() {
  std::vector<std::string> seed_table;
  // digits.txt is the the RAND table of random digits available from
  // http://www.rand.org/content/dam/rand/pubs/monograph_reports/MR1418/MR1418.digits.txt.zip
  std::ifstream seed_file("digits.txt");
  if (seed_file.good()) {
    std::string s;
    while (std::getline(seed_file, s)) {
      std::istringstream is(s);
      is >> s;
      std::string seed;
      for (int i = 0; i < 10; ++i) {
        is >> s;
        seed += s;
      }
      seed_table.push_back(seed);
    }
  } else {
    // lines 9077 - 9086, in case the full table is not available
    seed_table.push_back("91486866851718646640939841418208696361599740927130");
    seed_table.push_back("27085545979733960919843574531542387910433856492503");
    seed_table.push_back("50144629743871836759770203434855225597834197811609");
    seed_table.push_back("06513031977786096289011670052044879745609312522650");
    seed_table.push_back("27360659086699793904300006309496840459265099570954");
    seed_table.push_back("02069431357777418334914891893312167990996169079700");
    seed_table.push_back("32331134586674046873033847880723099262580390063311");
    seed_table.push_back("94820103089730893927855765220910774260125296553336");
    seed_table.push_back("13206851466217116627804451376819961678521211905382");
    seed_table.push_back("73582504716450282034788914850174375938706490445325");
  }
  exrandom::rand_table D;
  exrandom::u_rand<exrandom::rand_table> x(D);
  exrandom::unit_normal_dist<exrandom::rand_table> U(D);
  exrandom::table_gen g;
  for (size_t k = 0; k < seed_table.size(); ++k) {
    try {
      g.seed(seed_table[k]);
      long long c0 = D.count();
      U.generate(g,x);
      long long c1 = D.count();
      std::string x0 = x.print();
      std::string x1 = x.print_fixed(g, 6);
      std::string x2 = x.print();
      long long c2 = D.count();
      std::cout << seed_table[k].substr(0,int(c1-c0)) << "|"
                << seed_table[k].substr(int(c1-c0), int(c2-c1)) << "... -> ";
      std::cout << x0 << " = " << x2 << " ~ " << x1 << "\n";
    }
    catch (const std::exception&) {
      std::cerr << seed_table[k] << "... overflow\n";
    }
  }
}
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60

multiprec_test

#include 
#include 
#if !defined(EXRANDOM_HAVE_MPREAL)
#define EXRANDOM_HAVE_MPREAL 0
#endif
#if !defined(EXRANDOM_HAVE_BOOST_DEC)
#define EXRANDOM_HAVE_BOOST_DEC 0
#endif
#if !defined(EXRANDOM_HAVE_BOOST_GMP)
#define EXRANDOM_HAVE_BOOST_GMP 0
#endif
#if !defined(EXRANDOM_HAVE_BOOST_MPFR)
#define EXRANDOM_HAVE_BOOST_MPFR 0
#endif
#if !defined(EXRANDOM_HAVE_BOOST_FLOAT128)
#define EXRANDOM_HAVE_BOOST_FLOAT128 0
#endif
#if EXRANDOM_HAVE_MPREAL
#include 
#endif
#if EXRANDOM_HAVE_BOOST_MPFR
#include 
#endif
#if EXRANDOM_HAVE_BOOST_GMP
#include 
#endif
#if EXRANDOM_HAVE_BOOST_DEC
#include 
#endif
#if EXRANDOM_HAVE_BOOST_FLOAT128
#include 
#endif
// N.B. This #include should come *after* the headers which define the
// floating point types.
#include 
int main() {
#if EXRANDOM_HAVE_MPREAL || EXRANDOM_HAVE_BOOST_MPFR || \
  EXRANDOM_HAVE_BOOST_GMP || EXRANDOM_HAVE_BOOST_DEC || \
  EXRANDOM_HAVE_BOOST_FLOAT128
  unsigned s = std::random_device()(); // Set seed from random_device
  std::mt19937 g(s);                   // Initialize URNG
  std::cout << "Seed set to " << s << "\n\n";
#endif
#if EXRANDOM_HAVE_MPREAL
  {
    typedef mpfr::mpreal real;
    real::set_default_prec(240);
    exrandom::unit_normal_distribution<real> N;
    std::cout
      << "Sampling mpreal from the unit normal distribution:\n"
      << std::setprecision(mpfr::bits2digits(real::get_default_prec()));
    for (int i = 0; i < 20; ++i)
      std::cout << N(g) << "\n";
    std::cout << "\n";
  }
#else
  std::cerr << "This example requires mpreal.h\n\n";
#endif
#if EXRANDOM_HAVE_BOOST_MPFR
  {
    typedef boost::multiprecision::mpfr_float_50 real;
    exrandom::unit_normal_distribution<real> N;
    std::cout
      << "Sampling mpfr_float_50 from the unit normal distribution:\n"
      << std::setprecision(50);
    for (int i = 0; i < 20; ++i)
      std::cout << N(g) << "\n";
    std::cout << "\n";
  }
#else
  std::cerr << "This example requires boost/multiprecision/mpfr.hpp\n\n";
#endif
#if EXRANDOM_HAVE_BOOST_GMP
  {
    typedef boost::multiprecision::mpf_float_50 real;
    exrandom::unit_normal_distribution<real> N;
    std::cout
      << "Sampling mpf_float_50 from the unit normal distribution:\n"
      << std::setprecision(50);
    for (int i = 0; i < 20; ++i)
      std::cout << N(g) << "\n";
    std::cout << "\n";
  }
#else
  std::cerr << "This example requires boost/multiprecision/gmp.hpp\n\n";
#endif
#if EXRANDOM_HAVE_BOOST_DEC
  {
    typedef
      boost::multiprecision::number<boost::multiprecision::cpp_dec_float<20U>>
      real;
    exrandom::unit_normal_distribution<real> N;
    std::cout <<
      "Sampling number> from the unit normal distribution:\n"
      << std::setprecision(20);
    for (int i = 0; i < 20; ++i)
      std::cout << N(g) << "\n";
    std::cout << "\n";
  }
#else
  std::cerr
    << "This example requires boost/multiprecision/cpp_dec_float.hpp\n\n";
#endif
#if EXRANDOM_HAVE_BOOST_FLOAT128
  {
    typedef boost::multiprecision::float128 real;
    exrandom::unit_normal_distribution<real> N;
    std::cout
      << "Sampling float128 exactly from the unit normal distribution:\n"
      << std::setprecision(33);
    for (int i = 0; i < 20; ++i)
      std::cout << N(g) << "\n";
    std::cout << "\n";
  }
#else
  std::cerr << "This example requires boost/multiprecision/float128.hpp\n\n";
#endif
}
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118

hist_data

#include 
#include 
#include 
#include 
#include 
#include 
#include 
int main() {
  unsigned s = std::random_device()(); // Set seed from random_device
  std::mt19937 g(s);                   // Initialize URNG
  std::cerr << "Seed set to " << s << "\n"
            << "columns are l n prob(n/2^l)\n\n";
  const unsigned b = 2;
  exrandom::rand_digit<b> D;
  exrandom::unit_normal_dist<exrandom::rand_digit<b>> N(D);
  exrandom::u_rand<exrandom::rand_digit<b>> x(N.digit_generator());
  long long num = 5000000LL;
  unsigned maxl = 8, maxn = 3;
  std::map<std::pair<long long, long long>, long long> hist;
  for (long long i = 0; i < num; ++i) {
    N.generate(g, x);
    if (!(x.ndigits() <= maxl && x.integer() < maxn))
      continue;
    if (x.sign() < 1) x.negate();
    ++hist[x.rational()];
  }
  std::cout << std::setprecision(10);
  for (long l = 0, d = 1; l <= long(maxl); ++l, d *= 2)
    for (long n = 0; n < long(maxn)*d; ++n)
      std::cout << l << " " << n << " "
                << hist[std::pair<long long,long long>(n,d)]/double(num)
                << "\n";
}
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33

Seed set to 2659261496
columns are l n prob(n/2^l)

0 0 0.398886
0 1 0
0 2 0
1 0 0.0995716
1 1 0.0028946
1 2 0.0907916
1 3 0.013058
1 4 0.0101286
1 5 0.0003244
2 0 0.0250126
2 1 0.0125894
2 2 0.0374668
2 3 0.0054856
2 4 0.0339686
2 5 0.0066612
2 6 0.023425
2 7 0.0078724
2 8 0.0082116
2 9 0.0010498
2 10 0.0022626
2 11 0.0005622
3 0 0.0062726
3 1 0.0046628
3 2 0.007737
3 3 0.0032878
3 4 0.0093202
3 5 0.002567
3 6 0.011083
3 7 0.0028742
3 8 0.009916
3 9 0.0025964
3 10 0.0092066
3 11 0.002766
3 12 0.0080532
3 13 0.0028834
3 14 0.0063594
3 15 0.0027322
3 16 0.003031
3 17 0.0004968
3 18 0.002043
3 19 0.000556
3 20 0.001284
3 21 0.0004462
3 22 0.0007146
3 23 0.0002926
4 0 0.0015626
4 1 0.001343
4 2 0.0017742
4 3 0.0011566
4 4 0.0019266
4 5 0.001005
4 6 0.0021374
4 7 0.000878
4 8 0.002337
4 9 0.0008076
4 10 0.0025204
4 11 0.0007606
4 12 0.002736
4 13 0.0008342
4 14 0.0029704
4 15 0.0009782
4 16 0.0027164
4 17 0.0007852
4 18 0.0026088
4 19 0.000789
4 20 0.0025176
4 21 0.0008134
4 22 0.002436
4 23 0.000829
4 24 0.0022638
4 25 0.0008578
4 26 0.0021186
4 27 0.0008554
4 28 0.0018986
4 29 0.0008566
4 30 0.001706
4 31 0.0008322
4 32 0.0009082
4 33 0.0001574
4 34 0.0007646
4 35 0.0001898
4 36 0.000643
4 37 0.0001792
4 38 0.0005344
4 39 0.0001746
4 40 0.000441
4 41 0.0001678
4 42 0.0003568
4 43 0.0001572
4 44 0.0002692
4 45 0.000119
4 46 0.0002048
4 47 9.88e-05
5 0 0.0003854
5 1 0.0003664
5 2 0.0004178
5 3 0.0003428
5 4 0.0004458
5 5 0.000311
5 6 0.0004526
5 7 0.0002892
5 8 0.0005004
5 9 0.0002944
5 10 0.0005306
5 11 0.0002478
5 12 0.0005344
5 13 0.000251
5 14 0.0005522
5 15 0.0002252
5 16 0.0005806
5 17 0.0002168
5 18 0.0006126
5 19 0.0002154
5 20 0.0006336
5 21 0.0002118
5 22 0.0006498
5 23 0.0002134
5 24 0.000698
5 25 0.0002324
5 26 0.0007056
5 27 0.000246
5 28 0.0007448
5 29 0.0002612
5 30 0.0007674
5 31 0.000278
5 32 0.000667
5 33 0.0002182
5 34 0.000684
5 35 0.00022
5 36 0.0006806
5 37 0.0002102
5 38 0.0006718
5 39 0.0002118
5 40 0.000639
5 41 0.0002174
5 42 0.000641
5 43 0.0002288
5 44 0.000625
5 45 0.0002278
5 46 0.0006266
5 47 0.0002372
5 48 0.0006008
5 49 0.0002328
5 50 0.0005752
5 51 0.0002412
5 52 0.000552
5 53 0.0002436
5 54 0.0005406
5 55 0.000233
5 56 0.0005074
5 57 0.00024
5 58 0.0004752
5 59 0.0002396
5 60 0.0004838
5 61 0.000216
5 62 0.0004584
5 63 0.0002322
5 64 0.0002296
5 65 4.94e-05
5 66 0.0002218
5 67 5.38e-05
5 68 0.0002056
5 69 5.06e-05
5 70 0.000209
5 71 5.66e-05
5 72 0.0001856
5 73 5.04e-05
5 74 0.0001632
5 75 4.62e-05
5 76 0.0001462
5 77 4.78e-05
5 78 0.0001416
5 79 4.82e-05
5 80 0.000125
5 81 4.62e-05
5 82 0.000115
5 83 4.42e-05
5 84 9.42e-05
5 85 4.74e-05
5 86 9.2e-05
5 87 3.72e-05
5 88 8.3e-05
5 89 4.3e-05
5 90 6.88e-05
5 91 3.98e-05
5 92 6.62e-05
5 93 4.02e-05
5 94 5.14e-05
5 95 2.96e-05
6 0 9.58e-05
6 1 9.14e-05
6 2 8.9e-05
6 3 8.68e-05
6 4 0.000101
6 5 9.1e-05
6 6 0.0001034
6 7 8.58e-05
6 8 0.000107
6 9 8.52e-05
6 10 0.0001272
6 11 8.06e-05
6 12 0.0001208
6 13 8.38e-05
6 14 0.0001262
6 15 7.34e-05
6 16 0.0001142
6 17 7.1e-05
6 18 0.0001272
6 19 6.32e-05
6 20 0.000131
6 21 6.86e-05
6 22 0.0001356
6 23 6.86e-05
6 24 0.0001392
6 25 6.52e-05
6 26 0.0001418
6 27 7.16e-05
6 28 0.0001424
6 29 5.78e-05
6 30 0.0001368
6 31 6.34e-05
6 32 0.0001466
6 33 5.98e-05
6 34 0.0001484
6 35 5.92e-05
6 36 0.000153
6 37 5.5e-05
6 38 0.0001504
6 39 5.76e-05
6 40 0.0001644
6 41 5.52e-05
6 42 0.0001682
6 43 5.22e-05
6 44 0.0001722
6 45 6.06e-05
6 46 0.0001656
6 47 5.7e-05
6 48 0.0001786
6 49 5.62e-05
6 50 0.0001736
6 51 5.92e-05
6 52 0.0001768
6 53 5.78e-05
6 54 0.0001746
6 55 6.38e-05
6 56 0.0001902
6 57 6.12e-05
6 58 0.000191
6 59 6.32e-05
6 60 0.000206
6 61 7.14e-05
6 62 0.0001982
6 63 7.56e-05
6 64 0.0001836
6 65 5.5e-05
6 66 0.0001696
6 67 5.4e-05
6 68 0.0001736
6 69 5.96e-05
6 70 0.000182
6 71 6.14e-05
6 72 0.0001786
6 73 5.1e-05
6 74 0.000157
6 75 5.12e-05
6 76 0.0001632
6 77 5.46e-05
6 78 0.0001668
6 79 5.48e-05
6 80 0.0001714
6 81 5.86e-05
6 82 0.00017
6 83 5.84e-05
6 84 0.0001594
6 85 6.38e-05
6 86 0.0001668
6 87 6.56e-05
6 88 0.0001586
6 89 6.36e-05
6 90 0.0001606
6 91 5.76e-05
6 92 0.0001566
6 93 5.6e-05
6 94 0.0001484
6 95 6.32e-05
6 96 0.0001644
6 97 5.66e-05
6 98 0.0001496
6 99 5.9e-05
6 100 0.000148
6 101 6.02e-05
6 102 0.0001458
6 103 6e-05
6 104 0.0001438
6 105 6.12e-05
6 106 0.0001456
6 107 6.5e-05
6 108 0.0001382
6 109 6.4e-05
6 110 0.0001322
6 111 6.36e-05
6 112 0.0001334
6 113 6.38e-05
6 114 0.0001366
6 115 5.92e-05
6 116 0.0001236
6 117 5.98e-05
6 118 0.0001248
6 119 5.98e-05
6 120 0.000128
6 121 5.72e-05
6 122 0.0001198
6 123 5.42e-05
6 124 0.0001158
6 125 6.12e-05
6 126 0.0001102
6 127 6.88e-05
6 128 5.98e-05
6 129 9.8e-06
6 130 6.44e-05
6 131 1.36e-05
6 132 5.3e-05
6 133 1.38e-05
6 134 5.88e-05
6 135 1.18e-05
6 136 5.28e-05
6 137 1.08e-05
6 138 5.52e-05
6 139 1.06e-05
6 140 4.96e-05
6 141 1.22e-05
6 142 4.52e-05
6 143 1.36e-05
6 144 4.76e-05
6 145 1.12e-05
6 146 4.26e-05
6 147 1.34e-05
6 148 4.66e-05
6 149 1.28e-05
6 150 4.02e-05
6 151 1.62e-05
6 152 4.2e-05
6 153 1.24e-05
6 154 3.68e-05
6 155 1.34e-05
6 156 3.56e-05
6 157 1.2e-05
6 158 3.82e-05
6 159 1.46e-05
6 160 3.62e-05
6 161 1.24e-05
6 162 3e-05
6 163 1.38e-05
6 164 2.98e-05
6 165 1.1e-05
6 166 3e-05
6 167 1.2e-05
6 168 3.2e-05
6 169 1.32e-05
6 170 2.8e-05
6 171 1.12e-05
6 172 2.54e-05
6 173 1.16e-05
6 174 2.04e-05
6 175 1.06e-05
6 176 2.22e-05
6 177 8.6e-06
6 178 1.62e-05
6 179 8e-06
6 180 1.9e-05
6 181 8.4e-06
6 182 2.14e-05
6 183 8.6e-06
6 184 1.64e-05
6 185 9.8e-06
6 186 1.46e-05
6 187 8.4e-06
6 188 1.94e-05
6 189 7.4e-06
6 190 1.58e-05
6 191 6.6e-06
7 0 2.16e-05
7 1 2.42e-05
7 2 2.12e-05
7 3 2.5e-05
7 4 2.4e-05
7 5 2.52e-05
7 6 2.96e-05
7 7 2.06e-05
7 8 2.72e-05
7 9 2.32e-05
7 10 2.66e-05
7 11 1.94e-05
7 12 3.26e-05
7 13 2.2e-05
7 14 2.54e-05
7 15 2.64e-05
7 16 2.86e-05
7 17 1.96e-05
7 18 3.26e-05
7 19 1.92e-05
7 20 3.04e-05
7 21 2.28e-05
7 22 2.8e-05
7 23 2.28e-05
7 24 2.98e-05
7 25 1.84e-05
7 26 3.58e-05
7 27 1.54e-05
7 28 2.86e-05
7 29 1.76e-05
7 30 2.96e-05
7 31 1.92e-05
7 32 2.56e-05
7 33 2.2e-05
7 34 2.94e-05
7 35 1.76e-05
7 36 2.74e-05
7 37 1.8e-05
7 38 3.6e-05
7 39 1.64e-05
7 40 2.94e-05
7 41 1.54e-05
7 42 3.38e-05
7 43 1.32e-05
7 44 3.72e-05
7 45 1.5e-05
7 46 3.12e-05
7 47 1.82e-05
7 48 3.02e-05
7 49 1.82e-05
7 50 3.5e-05
7 51 1.26e-05
7 52 3.18e-05
7 53 1.26e-05
7 54 3.64e-05
7 55 1.38e-05
7 56 3.58e-05
7 57 1.48e-05
7 58 2.96e-05
7 59 1.62e-05
7 60 3.68e-05
7 61 1.5e-05
7 62 3.5e-05
7 63 1.48e-05
7 64 3.34e-05
7 65 1.68e-05
7 66 3.58e-05
7 67 1.32e-05
7 68 4.12e-05
7 69 1.38e-05
7 70 3.56e-05
7 71 1.44e-05
7 72 3.98e-05
7 73 1.22e-05
7 74 3.94e-05
7 75 1.26e-05
7 76 3.62e-05
7 77 1.68e-05
7 78 4.02e-05
7 79 1.38e-05
7 80 3.98e-05
7 81 1.78e-05
7 82 3.84e-05
7 83 1.18e-05
7 84 4.22e-05
7 85 1.42e-05
7 86 3.82e-05
7 87 1.06e-05
7 88 4.62e-05
7 89 1.74e-05
7 90 3.86e-05
7 91 1.24e-05
7 92 4e-05
7 93 1.68e-05
7 94 4.6e-05
7 95 1.52e-05
7 96 4.56e-05
7 97 1.52e-05
7 98 4.48e-05
7 99 1.74e-05
7 100 4.34e-05
7 101 1.62e-05
7 102 4.36e-05
7 103 1.88e-05
7 104 4.4e-05
7 105 1.42e-05
7 106 4.58e-05
7 107 1.96e-05
7 108 4.56e-05
7 109 1.46e-05
7 110 4.94e-05
7 111 1.74e-05
7 112 4.58e-05
7 113 1.6e-05
7 114 4.92e-05
7 115 1.64e-05
7 116 4.48e-05
7 117 1.84e-05
7 118 5.18e-05
7 119 1.82e-05
7 120 4.6e-05
7 121 1.8e-05
7 122 5.3e-05
7 123 2.22e-05
7 124 5.22e-05
7 125 1.98e-05
7 126 4.78e-05
7 127 1.92e-05
7 128 4.62e-05
7 129 1.5e-05
7 130 3.84e-05
7 131 1.78e-05
7 132 4.3e-05
7 133 1.34e-05
7 134 4.62e-05
7 135 1e-05
7 136 4.84e-05
7 137 1.04e-05
7 138 4.5e-05
7 139 1.1e-05
7 140 4.18e-05
7 141 1.58e-05
7 142 4.56e-05
7 143 1.4e-05
7 144 3.94e-05
7 145 1.46e-05
7 146 4.28e-05
7 147 1.76e-05
7 148 4.64e-05
7 149 1.58e-05
7 150 4.24e-05
7 151 1.78e-05
7 152 4.74e-05
7 153 1.46e-05
7 154 4.52e-05
7 155 1.34e-05
7 156 4.12e-05
7 157 1.6e-05
7 158 4.9e-05
7 159 1.76e-05
7 160 4.14e-05
7 161 1.44e-05
7 162 4.7e-05
7 163 1.48e-05
7 164 4.36e-05
7 165 1.56e-05
7 166 4.64e-05
7 167 1.62e-05
7 168 4.22e-05
7 169 1.44e-05
7 170 4.16e-05
7 171 1.86e-05
7 172 3.66e-05
7 173 1.28e-05
7 174 4.4e-05
7 175 1.46e-05
7 176 4.38e-05
7 177 1.3e-05
7 178 4.12e-05
7 179 1.5e-05
7 180 4.38e-05
7 181 1.3e-05
7 182 4.04e-05
7 183 1.5e-05
7 184 4e-05
7 185 1.6e-05
7 186 4.16e-05
7 187 1.38e-05
7 188 3.76e-05
7 189 1.5e-05
7 190 4.26e-05
7 191 1.5e-05
7 192 4.02e-05
7 193 1.82e-05
7 194 3.3e-05
7 195 1.8e-05
7 196 3.42e-05
7 197 1.7e-05
7 198 3.98e-05
7 199 1.56e-05
7 200 3.52e-05
7 201 1.56e-05
7 202 3.66e-05
7 203 1.56e-05
7 204 3.48e-05
7 205 1.68e-05
7 206 3.88e-05
7 207 1.74e-05
7 208 3.76e-05
7 209 1.44e-05
7 210 3.7e-05
7 211 1.38e-05
7 212 3.5e-05
7 213 1.6e-05
7 214 3.86e-05
7 215 1.68e-05
7 216 3.22e-05
7 217 1.86e-05
7 218 3.14e-05
7 219 1.98e-05
7 220 3.88e-05
7 221 1.66e-05
7 222 3.38e-05
7 223 1.46e-05
7 224 3.32e-05
7 225 1.3e-05
7 226 3.26e-05
7 227 1.52e-05
7 228 3.54e-05
7 229 1.66e-05
7 230 2.98e-05
7 231 1.7e-05
7 232 3.46e-05
7 233 2.1e-05
7 234 3.32e-05
7 235 1.34e-05
7 236 3.02e-05
7 237 1.5e-05
7 238 3.58e-05
7 239 1.62e-05
7 240 3.12e-05
7 241 1.42e-05
7 242 2.86e-05
7 243 1.66e-05
7 244 2.98e-05
7 245 1.8e-05
7 246 2.96e-05
7 247 1.62e-05
7 248 3.02e-05
7 249 1.34e-05
7 250 2.96e-05
7 251 1.94e-05
7 252 2.94e-05
7 253 1.56e-05
7 254 2.52e-05
7 255 1.46e-05
7 256 1.64e-05
7 257 3.2e-06
7 258 1.66e-05
7 259 2.6e-06
7 260 1.28e-05
7 261 3e-06
7 262 1.44e-05
7 263 5.2e-06
7 264 1.34e-05
7 265 3e-06
7 266 1.2e-05
7 267 5e-06
7 268 1.32e-05
7 269 3.8e-06
7 270 1.44e-05
7 271 4.8e-06
7 272 1.24e-05
7 273 4.6e-06
7 274 1.22e-05
7 275 4.4e-06
7 276 1.38e-05
7 277 3.6e-06
7 278 1.46e-05
7 279 4.4e-06
7 280 1.18e-05
7 281 3.2e-06
7 282 1e-05
7 283 3.6e-06
7 284 1.32e-05
7 285 3.8e-06
7 286 1.3e-05
7 287 3.6e-06
7 288 1.36e-05
7 289 3.2e-06
7 290 1.3e-05
7 291 4.2e-06
7 292 1.16e-05
7 293 3.8e-06
7 294 1.04e-05
7 295 3.4e-06
7 296 1.18e-05
7 297 4e-06
7 298 9.6e-06
7 299 3.4e-06
7 300 1.26e-05
7 301 3.8e-06
7 302 1.08e-05
7 303 4.2e-06
7 304 1.02e-05
7 305 5.6e-06
7 306 1.3e-05
7 307 4.2e-06
7 308 1.1e-05
7 309 3.4e-06
7 310 1.08e-05
7 311 3.6e-06
7 312 8e-06
7 313 2.8e-06
7 314 1.04e-05
7 315 2.4e-06
7 316 8.2e-06
7 317 4e-06
7 318 7.6e-06
7 319 4.4e-06
7 320 9.2e-06
7 321 3.4e-06
7 322 6.4e-06
7 323 3.2e-06
7 324 5e-06
7 325 2.6e-06
7 326 7.4e-06
7 327 2.6e-06
7 328 7.8e-06
7 329 3e-06
7 330 8.8e-06
7 331 4e-06
7 332 8.2e-06
7 333 3.8e-06
7 334 7e-06
7 335 3.8e-06
7 336 7e-06
7 337 4.8e-06
7 338 5.6e-06
7 339 2.6e-06
7 340 5.4e-06
7 341 2.2e-06
7 342 5.6e-06
7 343 2.4e-06
7 344 5.8e-06
7 345 4.2e-06
7 346 5.6e-06
7 347 3.8e-06
7 348 6.4e-06
7 349 1.4e-06
7 350 5.6e-06
7 351 1.6e-06
7 352 4.8e-06
7 353 2.2e-06
7 354 4.4e-06
7 355 2.2e-06
7 356 6e-06
7 357 2.6e-06
7 358 5.4e-06
7 359 1.2e-06
7 360 4.8e-06
7 361 4.4e-06
7 362 7.4e-06
7 363 3e-06
7 364 4e-06
7 365 2.6e-06
7 366 3e-06
7 367 2.8e-06
7 368 3.8e-06
7 369 2.4e-06
7 370 2.2e-06
7 371 2.4e-06
7 372 3.8e-06
7 373 2.6e-06
7 374 3.8e-06
7 375 8e-07
7 376 2.8e-06
7 377 2e-06
7 378 4.2e-06
7 379 3e-06
7 380 3.8e-06
7 381 1e-06
7 382 3.4e-06
7 383 2.6e-06
8 0 5.4e-06
8 1 9.2e-06
8 2 5.6e-06
8 3 5.4e-06
8 4 7.2e-06
8 5 6e-06
8 6 7.2e-06
8 7 8.4e-06
8 8 6.6e-06
8 9 4e-06
8 10 6.6e-06
8 11 7.2e-06
8 12 5.6e-06
8 13 6.8e-06
8 14 5.6e-06
8 15 6.4e-06
8 16 7e-06
8 17 8.4e-06
8 18 5.8e-06
8 19 5.2e-06
8 20 7.4e-06
8 21 8.6e-06
8 22 7e-06
8 23 4.8e-06
8 24 5.8e-06
8 25 7.4e-06
8 26 8e-06
8 27 5.4e-06
8 28 6.4e-06
8 29 4.2e-06
8 30 8e-06
8 31 4e-06
8 32 6.4e-06
8 33 5.2e-06
8 34 7.4e-06
8 35 3.4e-06
8 36 5.4e-06
8 37 6.6e-06
8 38 5e-06
8 39 6e-06
8 40 7e-06
8 41 4.6e-06
8 42 5e-06
8 43 3.8e-06
8 44 5.6e-06
8 45 4e-06
8 46 7e-06
8 47 4e-06
8 48 6.4e-06
8 49 5e-06
8 50 8.2e-06
8 51 5.2e-06
8 52 8.4e-06
8 53 5.2e-06
8 54 7e-06
8 55 4.2e-06
8 56 6.6e-06
8 57 4.2e-06
8 58 4.4e-06
8 59 4.4e-06
8 60 6.8e-06
8 61 3.4e-06
8 62 6e-06
8 63 5.6e-06
8 64 8.4e-06
8 65 4.4e-06
8 66 8.2e-06
8 67 4.2e-06
8 68 5.4e-06
8 69 3.8e-06
8 70 7.4e-06
8 71 5.6e-06
8 72 8.2e-06
8 73 4.2e-06
8 74 8e-06
8 75 3.6e-06
8 76 9.6e-06
8 77 5e-06
8 78 5.6e-06
8 79 5.4e-06
8 80 9.2e-06
8 81 3.6e-06
8 82 8.2e-06
8 83 5.6e-06
8 84 8e-06
8 85 4.2e-06
8 86 9e-06
8 87 5.6e-06
8 88 7.8e-06
8 89 3.2e-06
8 90 7.6e-06
8 91 3.6e-06
8 92 1.04e-05
8 93 5.2e-06
8 94 8e-06
8 95 4.2e-06
8 96 9.8e-06
8 97 2.8e-06
8 98 7.8e-06
8 99 3.2e-06
8 100 8.2e-06
8 101 5e-06
8 102 7.4e-06
8 103 4e-06
8 104 9.2e-06
8 105 4.2e-06
8 106 1.02e-05
8 107 4.6e-06
8 108 8.4e-06
8 109 3.8e-06
8 110 7e-06
8 111 3.8e-06
8 112 9e-06
8 113 3.4e-06
8 114 9.2e-06
8 115 3.8e-06
8 116 9.6e-06
8 117 2.4e-06
8 118 1.1e-05
8 119 2.2e-06
8 120 9e-06
8 121 4.8e-06
8 122 7.8e-06
8 123 4.4e-06
8 124 6.6e-06
8 125 3.2e-06
8 126 1.08e-05
8 127 3.6e-06
8 128 8.8e-06
8 129 4e-06
8 130 8.2e-06
8 131 3.2e-06
8 132 1.02e-05
8 133 3.2e-06
8 134 8.8e-06
8 135 3.6e-06
8 136 8.4e-06
8 137 3e-06
8 138 7.2e-06
8 139 3.8e-06
8 140 8.4e-06
8 141 2.6e-06
8 142 1.1e-05
8 143 3.8e-06
8 144 9.2e-06
8 145 3.4e-06
8 146 8.4e-06
8 147 2.2e-06
8 148 1.06e-05
8 149 3e-06
8 150 8.2e-06
8 151 3.8e-06
8 152 1e-05
8 153 4e-06
8 154 1.16e-05
8 155 4e-06
8 156 1.18e-05
8 157 3.4e-06
8 158 1.22e-05
8 159 3.8e-06
8 160 9.4e-06
8 161 4.4e-06
8 162 1.16e-05
8 163 4e-06
8 164 8e-06
8 165 4.2e-06
8 166 9.6e-06
8 167 3.8e-06
8 168 1.12e-05
8 169 4.2e-06
8 170 1.04e-05
8 171 2.8e-06
8 172 8.6e-06
8 173 3.2e-06
8 174 1.18e-05
8 175 1.8e-06
8 176 8.8e-06
8 177 5.4e-06
8 178 9.8e-06
8 179 1.6e-06
8 180 1.04e-05
8 181 4.2e-06
8 182 9.2e-06
8 183 3.4e-06
8 184 1.16e-05
8 185 4.6e-06
8 186 9.6e-06
8 187 2.4e-06
8 188 1.06e-05
8 189 3.6e-06
8 190 1.08e-05
8 191 2.6e-06
8 192 9.4e-06
8 193 2.6e-06
8 194 1e-05
8 195 4.6e-06
8 196 1.12e-05
8 197 4.8e-06
8 198 1.14e-05
8 199 3.4e-06
8 200 1.24e-05
8 201 4.2e-06
8 202 1.12e-05
8 203 3.4e-06
8 204 1.18e-05
8 205 3.6e-06
8 206 1.22e-05
8 207 3.6e-06
8 208 1.02e-05
8 209 2e-06
8 210 1e-05
8 211 5.6e-06
8 212 9.2e-06
8 213 3.6e-06
8 214 1.18e-05
8 215 3.6e-06
8 216 9.4e-06
8 217 3.6e-06
8 218 8.8e-06
8 219 3.2e-06
8 220 1.1e-05
8 221 2.4e-06
8 222 1.4e-05
8 223 5.2e-06
8 224 1.28e-05
8 225 5.6e-06
8 226 1.06e-05
8 227 5.8e-06
8 228 1.14e-05
8 229 5.2e-06
8 230 1.28e-05
8 231 5e-06
8 232 1.44e-05
8 233 5e-06
8 234 1.4e-05
8 235 4.6e-06
8 236 1.24e-05
8 237 3.6e-06
8 238 1.16e-05
8 239 4.8e-06
8 240 1.02e-05
8 241 5.4e-06
8 242 1.7e-05
8 243 5.4e-06
8 244 1.14e-05
8 245 4.2e-06
8 246 8.6e-06
8 247 4.8e-06
8 248 1.06e-05
8 249 4.6e-06
8 250 1.32e-05
8 251 6.6e-06
8 252 1.24e-05
8 253 5.6e-06
8 254 1.2e-05
8 255 3e-06
8 256 9e-06
8 257 3.8e-06
8 258 1.06e-05
8 259 2.6e-06
8 260 1.08e-05
8 261 3e-06
8 262 9.8e-06
8 263 3.6e-06
8 264 1e-05
8 265 3.2e-06
8 266 1.08e-05
8 267 2.8e-06
8 268 1e-05
8 269 3.6e-06
8 270 9.8e-06
8 271 4.6e-06
8 272 1.12e-05
8 273 4.2e-06
8 274 1.16e-05
8 275 4e-06
8 276 1.22e-05
8 277 3e-06
8 278 9e-06
8 279 4e-06
8 280 1.26e-05
8 281 4.6e-06
8 282 1.16e-05
8 283 5.2e-06
8 284 1.1e-05
8 285 2.6e-06
8 286 1.06e-05
8 287 4e-06
8 288 1.12e-05
8 289 3e-06
8 290 9.8e-06
8 291 2.8e-06
8 292 1.16e-05
8 293 6.8e-06
8 294 1.08e-05
8 295 4e-06
8 296 1.14e-05
8 297 5.2e-06
8 298 1.24e-05
8 299 3.2e-06
8 300 1.12e-05
8 301 4.6e-06
8 302 9.8e-06
8 303 3.2e-06
8 304 1.08e-05
8 305 4e-06
8 306 9.6e-06
8 307 4e-06
8 308 1.24e-05
8 309 4.6e-06
8 310 1.08e-05
8 311 6e-06
8 312 1.06e-05
8 313 4e-06
8 314 8e-06
8 315 3.8e-06
8 316 1.08e-05
8 317 3.8e-06
8 318 8.8e-06
8 319 3.6e-06
8 320 1.54e-05
8 321 2.2e-06
8 322 9.4e-06
8 323 1.8e-06
8 324 1.12e-05
8 325 3.4e-06
8 326 9.6e-06
8 327 3e-06
8 328 9.4e-06
8 329 3.4e-06
8 330 1.1e-05
8 331 3.6e-06
8 332 9.8e-06
8 333 1.8e-06
8 334 1.14e-05
8 335 3e-06
8 336 8.4e-06
8 337 3.6e-06
8 338 1.04e-05
8 339 4e-06
8 340 9.2e-06
8 341 4e-06
8 342 9.2e-06
8 343 2.2e-06
8 344 9.8e-06
8 345 3.6e-06
8 346 1.06e-05
8 347 2.6e-06
8 348 1.18e-05
8 349 4e-06
8 350 1e-05
8 351 4.2e-06
8 352 1.36e-05
8 353 3.2e-06
8 354 1.2e-05
8 355 4.2e-06
8 356 1e-05
8 357 5e-06
8 358 1.16e-05
8 359 3.8e-06
8 360 8.6e-06
8 361 2.8e-06
8 362 1.08e-05
8 363 3.6e-06
8 364 9e-06
8 365 3.8e-06
8 366 1.02e-05
8 367 4.6e-06
8 368 9.6e-06
8 369 4e-06
8 370 8.8e-06
8 371 4.4e-06
8 372 9.2e-06
8 373 3.8e-06
8 374 9.8e-06
8 375 4.4e-06
8 376 1e-05
8 377 4e-06
8 378 1.1e-05
8 379 4.6e-06
8 380 9.2e-06
8 381 4.4e-06
8 382 8.6e-06
8 383 2.8e-06
8 384 1.08e-05
8 385 4.6e-06
8 386 1.06e-05
8 387 3.8e-06
8 388 1e-05
8 389 3.4e-06
8 390 1.22e-05
8 391 3.6e-06
8 392 1.14e-05
8 393 3.4e-06
8 394 9.4e-06
8 395 3.8e-06
8 396 1.22e-05
8 397 4.4e-06
8 398 9.8e-06
8 399 2.2e-06
8 400 8.2e-06
8 401 4.2e-06
8 402 9.4e-06
8 403 5.4e-06
8 404 1.1e-05
8 405 4e-06
8 406 1e-05
8 407 5e-06
8 408 1.08e-05
8 409 4.6e-06
8 410 8.6e-06
8 411 5e-06
8 412 7.6e-06
8 413 5.2e-06
8 414 1.18e-05
8 415 2.8e-06
8 416 8.8e-06
8 417 6.2e-06
8 418 8.6e-06
8 419 3.8e-06
8 420 7.6e-06
8 421 3.6e-06
8 422 9.2e-06
8 423 3.6e-06
8 424 9e-06
8 425 3.6e-06
8 426 1.04e-05
8 427 4.8e-06
8 428 1.02e-05
8 429 3.8e-06
8 430 8e-06
8 431 3.2e-06
8 432 9.8e-06
8 433 2.8e-06
8 434 8.8e-06
8 435 1.6e-06
8 436 8.8e-06
8 437 4.4e-06
8 438 9.6e-06
8 439 3.6e-06
8 440 7e-06
8 441 3e-06
8 442 8.2e-06
8 443 5.4e-06
8 444 7e-06
8 445 2.8e-06
8 446 7e-06
8 447 3.2e-06
8 448 9.8e-06
8 449 3.8e-06
8 450 8.4e-06
8 451 5.4e-06
8 452 9.2e-06
8 453 4.2e-06
8 454 8.6e-06
8 455 3.6e-06
8 456 7.4e-06
8 457 4.8e-06
8 458 7.4e-06
8 459 3.4e-06
8 460 7.8e-06
8 461 4.6e-06
8 462 6.8e-06
8 463 4.8e-06
8 464 9.2e-06
8 465 2.6e-06
8 466 7.4e-06
8 467 3.8e-06
8 468 7.8e-06
8 469 4.8e-06
8 470 6.8e-06
8 471 2.8e-06
8 472 1.08e-05
8 473 2.2e-06
8 474 8.8e-06
8 475 4e-06
8 476 7.4e-06
8 477 2.8e-06
8 478 8.8e-06
8 479 5.8e-06
8 480 7.4e-06
8 481 3.4e-06
8 482 7.2e-06
8 483 5.4e-06
8 484 1.12e-05
8 485 3.4e-06
8 486 5.6e-06
8 487 3.4e-06
8 488 1.04e-05
8 489 3.2e-06
8 490 9.8e-06
8 491 5e-06
8 492 6e-06
8 493 3.4e-06
8 494 9.4e-06
8 495 2.2e-06
8 496 1.02e-05
8 497 5e-06
8 498 8.2e-06
8 499 3.6e-06
8 500 7.4e-06
8 501 3e-06
8 502 8.8e-06
8 503 3.8e-06
8 504 7.2e-06
8 505 4.4e-06
8 506 6e-06
8 507 4e-06
8 508 6.2e-06
8 509 3e-06
8 510 7e-06
8 511 5.6e-06
8 512 3.6e-06
8 513 6e-07
8 514 3.4e-06
8 515 6e-07
8 516 4.2e-06
8 517 6e-07
8 518 4.6e-06
8 519 6e-07
8 520 2.6e-06
8 521 6e-07
8 522 4.6e-06
8 523 0
8 524 2.8e-06
8 525 2e-07
8 526 4.4e-06
8 527 1e-06
8 528 4e-06
8 529 1.4e-06
8 530 4.4e-06
8 531 8e-07
8 532 3.4e-06
8 533 1.2e-06
8 534 2.4e-06
8 535 2e-07
8 536 3.4e-06
8 537 1e-06
8 538 4e-06
8 539 1.4e-06
8 540 3.2e-06
8 541 1e-06
8 542 3.8e-06
8 543 1.4e-06
8 544 3.6e-06
8 545 4e-07
8 546 5.2e-06
8 547 6e-07
8 548 3.2e-06
8 549 1e-06
8 550 3.2e-06
8 551 1.8e-06
8 552 3e-06
8 553 6e-07
8 554 3e-06
8 555 8e-07
8 556 3.8e-06
8 557 6e-07
8 558 4.4e-06
8 559 1.6e-06
8 560 2.2e-06
8 561 8e-07
8 562 2.2e-06
8 563 1e-06
8 564 3.8e-06
8 565 2e-07
8 566 2.4e-06
8 567 6e-07
8 568 3.4e-06
8 569 8e-07
8 570 2.8e-06
8 571 1e-06
8 572 3.2e-06
8 573 1.2e-06
8 574 1.8e-06
8 575 1e-06
8 576 2.2e-06
8 577 1.2e-06
8 578 2.6e-06
8 579 6e-07
8 580 2.8e-06
8 581 6e-07
8 582 1.6e-06
8 583 1.6e-06
8 584 2.4e-06
8 585 1.2e-06
8 586 4.2e-06
8 587 2e-07
8 588 3.2e-06
8 589 6e-07
8 590 3.2e-06
8 591 1e-06
8 592 3.6e-06
8 593 1.2e-06
8 594 3.4e-06
8 595 4e-07
8 596 2.2e-06
8 597 0
8 598 2.8e-06
8 599 6e-07
8 600 3e-06
8 601 1.6e-06
8 602 3.2e-06
8 603 1e-06
8 604 3e-06
8 605 8e-07
8 606 1.6e-06
8 607 0
8 608 2.6e-06
8 609 1e-06
8 610 2.8e-06
8 611 1e-06
8 612 2.2e-06
8 613 4e-07
8 614 2.6e-06
8 615 6e-07
8 616 1.8e-06
8 617 1.4e-06
8 618 2.2e-06
8 619 1.4e-06
8 620 2e-06
8 621 1.6e-06
8 622 2e-06
8 623 6e-07
8 624 1.8e-06
8 625 0
8 626 3.2e-06
8 627 1.2e-06
8 628 2.4e-06
8 629 6e-07
8 630 2e-06
8 631 6e-07
8 632 1e-06
8 633 4e-07
8 634 3e-06
8 635 1.8e-06
8 636 1.8e-06
8 637 1e-06
8 638 1.4e-06
8 639 4e-07
8 640 1.4e-06
8 641 4e-07
8 642 3.2e-06
8 643 8e-07
8 644 2e-06
8 645 1.2e-06
8 646 2.2e-06
8 647 4e-07
8 648 2.2e-06
8 649 4e-07
8 650 1.8e-06
8 651 8e-07
8 652 1.4e-06
8 653 1.6e-06
8 654 1.8e-06
8 655 6e-07
8 656 2.4e-06
8 657 1.4e-06
8 658 2e-06
8 659 6e-07
8 660 1.4e-06
8 661 6e-07
8 662 1.8e-06
8 663 1.4e-06
8 664 1.6e-06
8 665 1.4e-06
8 666 1.2e-06
8 667 6e-07
8 668 3e-06
8 669 8e-07
8 670 1.4e-06
8 671 0
8 672 8e-07
8 673 6e-07
8 674 1e-06
8 675 4e-07
8 676 1.8e-06
8 677 8e-07
8 678 1.6e-06
8 679 1.6e-06
8 680 2.6e-06
8 681 1e-06
8 682 2e-06
8 683 2e-07
8 684 2.4e-06
8 685 6e-07
8 686 2.4e-06
8 687 8e-07
8 688 1.2e-06
8 689 2e-07
8 690 1.6e-06
8 691 8e-07
8 692 1.6e-06
8 693 8e-07
8 694 1.2e-06
8 695 4e-07
8 696 6e-07
8 697 4e-07
8 698 1.8e-06
8 699 1.4e-06
8 700 1.8e-06
8 701 1.6e-06
8 702 1.2e-06
8 703 2e-07
8 704 1.4e-06
8 705 8e-07
8 706 6e-07
8 707 6e-07
8 708 1.8e-06
8 709 1.4e-06
8 710 1.4e-06
8 711 1e-06
8 712 2.2e-06
8 713 1.2e-06
8 714 8e-07
8 715 4e-07
8 716 1.2e-06
8 717 6e-07
8 718 1.2e-06
8 719 2e-07
8 720 8e-07
8 721 1e-06
8 722 1.2e-06
8 723 6e-07
8 724 1.2e-06
8 725 4e-07
8 726 8e-07
8 727 1.2e-06
8 728 1e-06
8 729 6e-07
8 730 1.4e-06
8 731 1.2e-06
8 732 1.8e-06
8 733 4e-07
8 734 8e-07
8 735 8e-07
8 736 1.6e-06
8 737 1.2e-06
8 738 2e-06
8 739 8e-07
8 740 6e-07
8 741 8e-07
8 742 1e-06
8 743 6e-07
8 744 1e-06
8 745 2e-07
8 746 1.2e-06
8 747 4e-07
8 748 4e-07
8 749 1.8e-06
8 750 1.8e-06
8 751 0
8 752 1.8e-06
8 753 4e-07
8 754 6e-07
8 755 2e-07
8 756 6e-07
8 757 2e-07
8 758 8e-07
8 759 4e-07
8 760 2e-07
8 761 8e-07
8 762 6e-07
8 763 6e-07
8 764 1.6e-06
8 765 0
8 766 1.6e-06
8 767 1.4e-06

discrete_count_bits

#include 
#include 
#include 
#include 
#include 
#include 
template<typename Generator>
void count_bits(Generator& g, long long num,
                int mu_num, int mu_den, int sigma_num, int sigma_den) {
  const unsigned b = 2;
  exrandom::rand_digit<b> D;
  exrandom::discrete_normal_dist<exrandom::rand_digit<b>>
    N(D, mu_num, mu_den, sigma_num, sigma_den);
  exrandom::i_rand<exrandom::rand_digit<b>> h(N.digit_generator());
  long long counta = 0, countb = 0, fractotal = 0;
  for (long long i = 0; i < num; ++i) {
    long long c0 = D.count();
    N.generate(g, h);
    counta += D.count() - c0;
    fractotal += h.entropy();
    h(g);
    countb += D.count() - c0;
  }
  double H;
  exrandom::aux_info::discrete_normal_norm(N.param(), H);
  H /= std::log(2.0);
  double B = counta/double(num),
    L = fractotal/double(num),
    F = countb/double(num),
    T = F - H;
  std::cout << "mu = " << N.mu_num() << "/" << N.mu_den()
            << ", sigma = " << N.sigma_num() << "/" << N.sigma_den()
            << ":\n";
  std::cout << "  = " << B
            << ",  = " << L
            << ",  =  +  = " << F
            << ",\n H = " << H
            << ", T =  - H = " << T << "" << std::endl;
}
int main() {
  unsigned s = std::random_device()(); // Set seed from random_device
  std::mt19937 g(s);                   // Initialize URNG
  std::cout << "Seed set to " << s << "\n\n";
  long long num = 2000000LL;
  std::cout
    << "Sampling from the discrete normal distribution (Algorithm D)\n";
  count_bits(g, num,  1, 7, 16, 100);
  count_bits(g, num,  2, 7, 16, 10);
  count_bits(g, num,  3, 7, 16, 1);
  count_bits(g, num,  0, 7, 160, 1);
  count_bits(g, num, -1, 7, 1600, 1);
  count_bits(g, num, -2, 7, 16000, 1);
  count_bits(g, num, -3, 7, 160000, 1);
  count_bits(g, num,  1, 7, 1600000, 1);
  count_bits(g, num,  0, 7, (1<<17)-1, 1);
  count_bits(g, num,  0, 7, (1<<17),   1);
  count_bits(g, num,  0, 1, (1<<17)+1, 1);
  count_bits(g, num,  0, 1, (1<<18)-1, 1);
  count_bits(g, num,  0, 1, (1<<18),   1);
  count_bits(g, num,  0, 1, (1<<18)+1, 1);
}
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61

#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
void print_hist(const std::map<int, long long>& hist) {
  const int maxlength = 70;     // Widest bar in histogram
  long long m = 0;
  for (auto p = hist.begin(); p != hist.end(); ++p)
    m = (std::max)(m, p->second);
  m = (std::max)(m / maxlength, 1LL);
  for (auto p = hist.begin(); p != hist.end(); ++p) {
    unsigned k = unsigned((p->second + m/2) / m);
    if (k) {
      std::cout << p->first << '\t';
      for (unsigned i = 0; i < k; ++i)
        std::cout << '*';
      std::cout << '\n';
    }
  }
  std::cout << std::endl;
}
int main() {
  unsigned s = std::random_device()(); // Set seed from random_device
  std::mt19937 g(s);                   // Initialize URNG
  std::cout << "Seed set to " << s << "\n" << std::endl;
  const unsigned b = 2;
  long long num = 1000000LL;
  std::cout << std::fixed << std::setprecision(5);
  {
    std::cout << "Sampling from unit normal distribution\n";
    exrandom::rand_digit<b> D;
    exrandom::unit_normal_dist<exrandom::rand_digit<b>> N(D);
    exrandom::u_rand<exrandom::rand_digit<b>>
      x(N.digit_generator());
    std::map<int, long long> bitcount;
    std::map<int, long long> fraccount;
    long long counta = 0, countb = 0, fractotal = 0;
    for (long long i = 0; i < num; ++i) {
      long long c0 = D.count();
      N.generate(g, x);
      counta += D.count() - c0;
      fractotal += x.ndigits();
      ++fraccount[int(x.ndigits())];
      x.value<double>(g);
      countb += D.count() - c0;
      ++bitcount[int(D.count() - c0)];
    }
    double H = exrandom::aux_info::normal_entropy() / std::log(2.0),
      Q = exrandom::aux_info::normal_mean_exponent(),
      B = counta/double(num),
      L = fractotal/double(num),
      C = B - L, T = C - H,
      F = countb/double(num);
    int d = std::numeric_limits<double>::digits + 1;
    std::cout << " = " << B
              << ",  = " << L
              << ", C =  -  = " << C
              << ",\nH = " << H
              << ", T = C - H = " << T
              << ",\nQ = " << Q
              << ", F = C - Q + " << d << " = " << C - Q + d
              << " (predicted) = " << F << " (measured)\n" << std::endl;
    std::cout << "Histogram of number of bits in fraction\n";
    print_hist(fraccount);
    std::cout << "Histogram of number of bits needed to generate doubles\n";
    print_hist(bitcount);
  }
  {
    std::cout << "Sampling from unit normal distribution (Kahn method)\n";
    exrandom::rand_digit<b> D;
    exrandom::unit_normal_kahn<exrandom::rand_digit<b>> N(D);
    exrandom::u_rand<exrandom::rand_digit<b>>
      x(N.digit_generator());
    std::map<int, long long> bitcount;
    std::map<int, long long> fraccount;
    long long counta = 0, countb = 0, fractotal = 0;
    for (long long i = 0; i < num; ++i) {
      long long c0 = D.count();
      N.generate(g, x);
      counta += D.count() - c0;
      fractotal += x.ndigits();
      ++fraccount[int(x.ndigits())];
      x.value<double>(g);
      countb += D.count() - c0;
      ++bitcount[int(D.count() - c0)];
    }
    double H = exrandom::aux_info::normal_entropy() / std::log(2.0),
      Q = exrandom::aux_info::normal_mean_exponent(),
      B = counta/double(num),
      L = fractotal/double(num),
      C = B - L, T = C - H,
      F = countb/double(num);
    int d = std::numeric_limits<double>::digits + 1;
    std::cout << " = " << B
              << ",  = " << L
              << ", C =  -  = " << C
              << ",\nH = " << H
              << ", T = C - H = " << T
              << ",\nQ = " << Q
              << ", F = C - Q + " << d << " = " << C - Q + d
              << " (predicted) = " << F << " (measured)\n" << std::endl;
    std::cout << "Histogram of number of bits in fraction\n";
    print_hist(fraccount);
    std::cout << "Histogram of number of bits needed to generate doubles\n";
    print_hist(bitcount);
  }
  {
    std::cout
      << "Sampling from unit exponential distribution (Algorithm V)\n";
    exrandom::rand_digit<b> D;
    exrandom::unit_exponential_dist<exrandom::rand_digit<b>,false> N(D);
    exrandom::u_rand<exrandom::rand_digit<b>>
      x(N.digit_generator());
    std::map<int, long long> bitcount;
    std::map<int, long long> fraccount;
    long long counta = 0, countb = 0, fractotal = 0;
    for (long long i = 0; i < num; ++i) {
      long long c0 = D.count();
      N.generate(g, x);
      counta += D.count() - c0;
      fractotal += x.ndigits();
      ++fraccount[int(x.ndigits())];
      x.value<double>(g);
      countb += D.count() - c0;
      ++bitcount[int(D.count() - c0)];
    }
    double H = exrandom::aux_info::exponential_entropy() / std::log(2.0),
      Q = exrandom::aux_info::exponential_mean_exponent(),
      B = counta/double(num),
      L = fractotal/double(num),
      C = B - L, T = C - H,
      F = countb/double(num);
    int d = std::numeric_limits<double>::digits + 1;
    std::cout << " = " << B
              << ",  = " << L
              << ", C =  -  = " << C
              << ",\nH = " << H
              << ", T = C - H = " << T
              << ",\nQ = " << Q
              << ", F = C - Q + " << d << " = " << C - Q + d
              << " (predicted) = " << F << " (measured)\n" << std::endl;
  }
  {
    std::cout
      << "Sampling from unit exponential distribution (Algorithm E)\n";
    exrandom::rand_digit<b> D;
    exrandom::unit_exponential_dist<exrandom::rand_digit<b>> N(D);
    exrandom::u_rand<exrandom::rand_digit<b>>
      x(N.digit_generator());
    std::map<int, long long> bitcount;
    std::map<int, long long> fraccount;
    long long counta = 0, countb = 0, fractotal = 0;
    for (long long i = 0; i < num; ++i) {
      long long c0 = D.count();
      N.generate(g, x);
      counta += D.count() - c0;
      fractotal += x.ndigits();
      ++fraccount[int(x.ndigits())];
      x.value<double>(g);
      countb += D.count() - c0;
      ++bitcount[int(D.count() - c0)];
    }
    double H = exrandom::aux_info::exponential_entropy() / std::log(2.0),
      Q = exrandom::aux_info::exponential_mean_exponent(),
      B = counta/double(num),
      L = fractotal/double(num),
      C = B - L, T = C - H,
      F = countb/double(num);
    int d = std::numeric_limits<double>::digits + 1;
    std::cout << " = " << B
              << ",  = " << L
              << ", C =  -  = " << C
              << ",\nH = " << H
              << ", T = C - H = " << T
              << ",\nQ = " << Q
              << ", F = C - Q + " << d << " = " << C - Q + d
              << " (predicted) = " << F << " (measured)\n" << std::endl;
    std::cout << "Histogram of number of bits in fraction\n";
    print_hist(fraccount);
    std::cout << "Histogram of number of bits needed to generate doubles\n";
    print_hist(bitcount);
  }
  {
    std::cout
      << "Sampling from unit uniform distribution\n";
    exrandom::rand_digit<b> D;
    exrandom::unit_uniform_dist<exrandom::rand_digit<b>> N(D);
    exrandom::u_rand<exrandom::rand_digit<b>> x(N.digit_generator());
    std::map<int, long long> bitcount;
    std::map<int, long long> fraccount;
    long long counta = 0, countb = 0, fractotal = 0;
    for (long long i = 0; i < num; ++i) {
      long long c0 = D.count();
      N.generate(g, x);
      counta += D.count() - c0;
      fractotal += x.ndigits();
      ++fraccount[int(x.ndigits())];
      x.value<double>(g);
      countb += D.count() - c0;
      ++bitcount[int(D.count() - c0)];
    }
    double H = exrandom::aux_info::uniform_entropy() / std::log(2.0),
      Q = exrandom::aux_info::uniform_mean_exponent(),
      B = counta/double(num),
      L = fractotal/double(num),
      C = B - L, T = C - H,
      F = countb/double(num);
    int d = std::numeric_limits<double>::digits + 1;
    std::cout << " = " << B
              << ",  = " << L
              << ", C =  -  = " << C
              << ",\nH = " << H
              << ", T = C - H = " << T
              << ",\nQ = " << Q
              << ", F = C - Q + " << d << " = " << C - Q + d
              << " (predicted) = " << F << " (measured)\n" << std::endl;
    std::cout << "Histogram of number of bits in fraction\n";
    print_hist(fraccount);
    std::cout << "Histogram of number of bits needed to generate doubles\n";
    print_hist(bitcount);
  }
}
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229