From: Paolo Carlini Date: Thu, 24 Mar 2011 16:49:19 +0000 (+0000) Subject: random.h (negative_binomial_distribution<>:: negative_binomial_distribution(_IntType... X-Git-Url: https://git.libre-soc.org/?a=commitdiff_plain;h=113b21bdaf45d2b80d1585b4b2c9cb20395667be;p=gcc.git random.h (negative_binomial_distribution<>:: negative_binomial_distribution(_IntType, double), [...]): Fix construction of _M_gd. 2011-03-24 Paolo Carlini * include/bits/random.h (negative_binomial_distribution<>:: negative_binomial_distribution(_IntType, double), negative_binomial_distribution<>:: negative_binomial_distribution(const param_type&)): Fix construction of _M_gd. * include/bits/random.tcc (negative_binomial_distribution<>:: operator()): Fix computation, per Leger's algorithm. * testsuite/util/testsuite_random.h (discrete_pdf, negative_binomial_pdf, poisson_pdf, uniform_int_pdf): New. (binomial_pdf): Swap last two parameters. * testsuite/26_numerics/random/discrete_distribution/ operators/values.cc: New. * testsuite/26_numerics/random/negative_binomial_distribution/ operators/values.cc: Likewise. * testsuite/26_numerics/random/poisson_distribution/ operators/values.cc: Likewise. * testsuite/26_numerics/random/uniform_int_distribution/ operators/values.cc: Likewise. * testsuite/26_numerics/random/binomial_distribution/ operators/values.cc: Adjust. From-SVN: r171411 --- diff --git a/libstdc++-v3/ChangeLog b/libstdc++-v3/ChangeLog index 538524880de..9182dc0ae4e 100644 --- a/libstdc++-v3/ChangeLog +++ b/libstdc++-v3/ChangeLog @@ -1,3 +1,26 @@ +2011-03-24 Paolo Carlini + + * include/bits/random.h (negative_binomial_distribution<>:: + negative_binomial_distribution(_IntType, double), + negative_binomial_distribution<>:: + negative_binomial_distribution(const param_type&)): Fix + construction of _M_gd. + * include/bits/random.tcc (negative_binomial_distribution<>:: + operator()): Fix computation, per Leger's algorithm. + * testsuite/util/testsuite_random.h (discrete_pdf, + negative_binomial_pdf, poisson_pdf, uniform_int_pdf): New. + (binomial_pdf): Swap last two parameters. + * testsuite/26_numerics/random/discrete_distribution/ + operators/values.cc: New. + * testsuite/26_numerics/random/negative_binomial_distribution/ + operators/values.cc: Likewise. + * testsuite/26_numerics/random/poisson_distribution/ + operators/values.cc: Likewise. + * testsuite/26_numerics/random/uniform_int_distribution/ + operators/values.cc: Likewise. + * testsuite/26_numerics/random/binomial_distribution/ + operators/values.cc: Adjust. + 2011-03-24 Rainer Orth * config/abi/post/solaris2.8/baseline_symbols.txt: Regenerate. diff --git a/libstdc++-v3/include/bits/random.h b/libstdc++-v3/include/bits/random.h index 26cec8a885e..8b09a98c37b 100644 --- a/libstdc++-v3/include/bits/random.h +++ b/libstdc++-v3/include/bits/random.h @@ -3611,8 +3611,7 @@ _GLIBCXX_BEGIN_NAMESPACE_VERSION param_type(double __p = 0.5) : _M_p(__p) { - _GLIBCXX_DEBUG_ASSERT((_M_p > 0.0) - && (_M_p < 1.0)); + _GLIBCXX_DEBUG_ASSERT((_M_p > 0.0) && (_M_p < 1.0)); _M_initialize(); } @@ -3782,7 +3781,9 @@ _GLIBCXX_BEGIN_NAMESPACE_VERSION explicit param_type(_IntType __k = 1, double __p = 0.5) : _M_k(__k), _M_p(__p) - { } + { + _GLIBCXX_DEBUG_ASSERT((_M_k > 0) && (_M_p > 0.0) && (_M_p <= 1.0)); + } _IntType k() const @@ -3803,12 +3804,12 @@ _GLIBCXX_BEGIN_NAMESPACE_VERSION explicit negative_binomial_distribution(_IntType __k = 1, double __p = 0.5) - : _M_param(__k, __p), _M_gd(__k, __p / (1.0 - __p)) + : _M_param(__k, __p), _M_gd(__k, 1.0) { } explicit negative_binomial_distribution(const param_type& __p) - : _M_param(__p), _M_gd(__p.k(), __p.p() / (1.0 - __p.p())) + : _M_param(__p), _M_gd(__p.k(), 1.0) { } /** diff --git a/libstdc++-v3/include/bits/random.tcc b/libstdc++-v3/include/bits/random.tcc index 4b17e914e56..0bbc9fd400d 100644 --- a/libstdc++-v3/include/bits/random.tcc +++ b/libstdc++-v3/include/bits/random.tcc @@ -1075,7 +1075,7 @@ _GLIBCXX_BEGIN_NAMESPACE_VERSION return __is; } - + // This is Leger's algorithm. template template typename negative_binomial_distribution<_IntType>::result_type @@ -1085,7 +1085,8 @@ _GLIBCXX_BEGIN_NAMESPACE_VERSION const double __y = _M_gd(__urng); // XXX Is the constructor too slow? - std::poisson_distribution __poisson(__y); + std::poisson_distribution __poisson(__y * (1.0 - p()) + / p()); return __poisson(__urng); } @@ -1099,10 +1100,10 @@ _GLIBCXX_BEGIN_NAMESPACE_VERSION typedef typename std::gamma_distribution::param_type param_type; - const double __y = - _M_gd(__urng, param_type(__p.k(), __p.p() / (1.0 - __p.p()))); + const double __y = _M_gd(__urng, param_type(__p.k(), 1.0)); - std::poisson_distribution __poisson(__y); + std::poisson_distribution __poisson(__y * (1.0 - __p.p()) + / __p.p() ); return __poisson(__urng); } diff --git a/libstdc++-v3/testsuite/26_numerics/random/binomial_distribution/operators/values.cc b/libstdc++-v3/testsuite/26_numerics/random/binomial_distribution/operators/values.cc index c0248ab914c..9fa8f033dff 100644 --- a/libstdc++-v3/testsuite/26_numerics/random/binomial_distribution/operators/values.cc +++ b/libstdc++-v3/testsuite/26_numerics/random/binomial_distribution/operators/values.cc @@ -33,16 +33,16 @@ void test01() std::binomial_distribution<> bd1(5, 0.3); auto bbd1 = std::bind(bd1, eng); - testDiscreteDist(bbd1, [](int n) { return binomial_pdf(n, 0.3, 5); } ); + testDiscreteDist(bbd1, [](int n) { return binomial_pdf(n, 5, 0.3); } ); std::binomial_distribution<> bd2(55, 0.3); auto bbd2 = std::bind(bd2, eng); - testDiscreteDist(bbd2, [](int n) { return binomial_pdf(n, 0.3, 55); } ); + testDiscreteDist(bbd2, [](int n) { return binomial_pdf(n, 55, 0.3); } ); // libstdc++/48114 std::binomial_distribution<> bd3(10, 0.75); auto bbd3 = std::bind(bd3, eng); - testDiscreteDist(bbd3, [](int n) { return binomial_pdf(n, 0.75, 10); } ); + testDiscreteDist(bbd3, [](int n) { return binomial_pdf(n, 10, 0.75); } ); } int main() diff --git a/libstdc++-v3/testsuite/26_numerics/random/discrete_distribution/operators/values.cc b/libstdc++-v3/testsuite/26_numerics/random/discrete_distribution/operators/values.cc new file mode 100644 index 00000000000..171087554da --- /dev/null +++ b/libstdc++-v3/testsuite/26_numerics/random/discrete_distribution/operators/values.cc @@ -0,0 +1,52 @@ +// { dg-options "-std=gnu++0x" } +// { dg-require-cstdint "" } +// +// Copyright (C) 2011 Free Software Foundation, Inc. +// +// This file is part of the GNU ISO C++ Library. This library is free +// software; you can redistribute it and/or modify it under the +// terms of the GNU General Public License as published by the +// Free Software Foundation; either version 3, or (at your option) +// any later version. +// +// This library is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +// GNU General Public License for more details. +// +// You should have received a copy of the GNU General Public License along +// with this library; see the file COPYING3. If not see +// . + +// 26.5.8.6.1 Class template discrete_distribution [rand.dist.samp.discrete] + +#include +#include +#include + +void test01() +{ + using namespace __gnu_test; + + std::mt19937 eng; + + std::discrete_distribution<> dd1({ }); + auto bdd1 = std::bind(dd1, eng); + testDiscreteDist(bdd1, [](int n) { return discrete_pdf(n, { }); } ); + + std::discrete_distribution<> dd2({ 1.0, 3.0, 2.0}); + auto bdd2 = std::bind(dd2, eng); + testDiscreteDist(bdd2, [](int n) + { return discrete_pdf(n, { 1.0, 3.0, 2.0}); } ); + + std::discrete_distribution<> dd3({ 2.0, 2.0, 1.0, 0.0, 4.0}); + auto bdd3 = std::bind(dd3, eng); + testDiscreteDist(bdd3, [](int n) + { return discrete_pdf(n, { 2.0, 2.0, 1.0, 0.0, 4.0}); } ); +} + +int main() +{ + test01(); + return 0; +} diff --git a/libstdc++-v3/testsuite/26_numerics/random/negative_binomial_distribution/operators/values.cc b/libstdc++-v3/testsuite/26_numerics/random/negative_binomial_distribution/operators/values.cc new file mode 100644 index 00000000000..fb5194aebca --- /dev/null +++ b/libstdc++-v3/testsuite/26_numerics/random/negative_binomial_distribution/operators/values.cc @@ -0,0 +1,55 @@ +// { dg-options "-std=gnu++0x" } +// { dg-require-cstdint "" } +// { dg-require-cmath "" } +// +// Copyright (C) 2011 Free Software Foundation, Inc. +// +// This file is part of the GNU ISO C++ Library. This library is free +// software; you can redistribute it and/or modify it under the +// terms of the GNU General Public License as published by the +// Free Software Foundation; either version 3, or (at your option) +// any later version. +// +// This library is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +// GNU General Public License for more details. +// +// You should have received a copy of the GNU General Public License along +// with this library; see the file COPYING3. If not see +// . + +// 26.5.8.3.4 Class template negative_binomial_distribution +// [rand.dist.bern.negbin] + +#include +#include +#include + +void test01() +{ + using namespace __gnu_test; + + std::mt19937 eng; + + std::negative_binomial_distribution<> nbd1(5, 0.3); + auto bnbd1 = std::bind(nbd1, eng); + testDiscreteDist(bnbd1, [](int n) + { return negative_binomial_pdf(n, 5, 0.3); } ); + + std::negative_binomial_distribution<> nbd2(55, 0.3); + auto bnbd2 = std::bind(nbd2, eng); + testDiscreteDist(bnbd2, [](int n) + { return negative_binomial_pdf(n, 55, 0.3); } ); + + std::negative_binomial_distribution<> nbd3(10, 0.75); + auto bnbd3 = std::bind(nbd3, eng); + testDiscreteDist(bnbd3, [](int n) + { return negative_binomial_pdf(n, 10, 0.75); } ); +} + +int main() +{ + test01(); + return 0; +} diff --git a/libstdc++-v3/testsuite/26_numerics/random/poisson_distribution/operators/values.cc b/libstdc++-v3/testsuite/26_numerics/random/poisson_distribution/operators/values.cc new file mode 100644 index 00000000000..adfe9b7528e --- /dev/null +++ b/libstdc++-v3/testsuite/26_numerics/random/poisson_distribution/operators/values.cc @@ -0,0 +1,51 @@ +// { dg-options "-std=gnu++0x" } +// { dg-require-cstdint "" } +// { dg-require-cmath "" } +// +// Copyright (C) 2011 Free Software Foundation, Inc. +// +// This file is part of the GNU ISO C++ Library. This library is free +// software; you can redistribute it and/or modify it under the +// terms of the GNU General Public License as published by the +// Free Software Foundation; either version 3, or (at your option) +// any later version. +// +// This library is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +// GNU General Public License for more details. +// +// You should have received a copy of the GNU General Public License along +// with this library; see the file COPYING3. If not see +// . + +// 26.5.8.4.1 Class template poisson_distribution [rand.dist.pois.poisson] + +#include +#include +#include + +void test01() +{ + using namespace __gnu_test; + + std::mt19937 eng; + + std::poisson_distribution<> pd1(3.0); + auto bpd1 = std::bind(pd1, eng); + testDiscreteDist(bpd1, [](int n) { return poisson_pdf(n, 3.0); } ); + + std::poisson_distribution<> pd2(15.0); + auto bpd2 = std::bind(pd2, eng); + testDiscreteDist(bpd2, [](int n) { return poisson_pdf(n, 15.0); } ); + + std::poisson_distribution<> pd3(30.0); + auto bpd3 = std::bind(pd3, eng); + testDiscreteDist(bpd3, [](int n) { return poisson_pdf(n, 30.0); } ); +} + +int main() +{ + test01(); + return 0; +} diff --git a/libstdc++-v3/testsuite/26_numerics/random/uniform_int_distribution/operators/values.cc b/libstdc++-v3/testsuite/26_numerics/random/uniform_int_distribution/operators/values.cc new file mode 100644 index 00000000000..8cb99a1a55d --- /dev/null +++ b/libstdc++-v3/testsuite/26_numerics/random/uniform_int_distribution/operators/values.cc @@ -0,0 +1,50 @@ +// { dg-options "-std=gnu++0x" } +// { dg-require-cstdint "" } +// +// Copyright (C) 2011 Free Software Foundation, Inc. +// +// This file is part of the GNU ISO C++ Library. This library is free +// software; you can redistribute it and/or modify it under the +// terms of the GNU General Public License as published by the +// Free Software Foundation; either version 3, or (at your option) +// any later version. +// +// This library is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +// GNU General Public License for more details. +// +// You should have received a copy of the GNU General Public License along +// with this library; see the file COPYING3. If not see +// . + +// 26.5.8.2.1 Class template uniform_int_distribution [rand.dist.uni.int] + +#include +#include +#include + +void test01() +{ + using namespace __gnu_test; + + std::mt19937 eng; + + std::uniform_int_distribution<> uid1(0, 2); + auto buid1 = std::bind(uid1, eng); + testDiscreteDist(buid1, [](int n) { return uniform_int_pdf(n, 0, 2); } ); + + std::uniform_int_distribution<> uid2(3, 7); + auto buid2 = std::bind(uid2, eng); + testDiscreteDist(buid2, [](int n) { return uniform_int_pdf(n, 3, 7); } ); + + std::uniform_int_distribution<> uid3(1, 20); + auto buid3 = std::bind(uid3, eng); + testDiscreteDist(buid3, [](int n) { return uniform_int_pdf(n, 1, 20); } ); +} + +int main() +{ + test01(); + return 0; +} diff --git a/libstdc++-v3/testsuite/util/testsuite_random.h b/libstdc++-v3/testsuite/util/testsuite_random.h index 417a95a8c6c..5b5575571d7 100644 --- a/libstdc++-v3/testsuite/util/testsuite_random.h +++ b/libstdc++-v3/testsuite/util/testsuite_random.h @@ -25,6 +25,7 @@ #define _GLIBCXX_TESTSUITE_RANDOM_H #include +#include #include namespace __gnu_test @@ -79,27 +80,27 @@ namespace __gnu_test else if (k == 1) return p; else - return 0; + return 0.0; } #ifdef _GLIBCXX_USE_C99_MATH_TR1 inline double - binomial_pdf(int k, double p, int n) + binomial_pdf(int k, int n, double p) { if (k < 0 || k > n) - return 0; + return 0.0; else { double q; - if (p == 0) - q = (k == 0) ? 1 : 0; - else if (p == 1) - q = (k == n) ? 1 : 0; + if (p == 0.0) + q = (k == 0) ? 1.0 : 0.0; + else if (p == 1.0) + q = (k == n) ? 1.0 : 0.0; else { - double ln_Cnk = (std::lgamma(n + 1) - std::lgamma(k + 1) - - std::lgamma(n - k + 1)); + double ln_Cnk = (std::lgamma(n + 1.0) - std::lgamma(k + 1.0) + - std::lgamma(n - k + 1.0)); q = ln_Cnk + k * std::log(p) + (n - k) * std::log1p(-p); q = std::exp(q); } @@ -109,16 +110,72 @@ namespace __gnu_test } #endif + inline double + discrete_pdf(int k, std::initializer_list wl) + { + if (!wl.size()) + wl = { 1.0 }; + + if (k < 0 || k >= wl.size()) + return 0.0; + else + { + double sum = 0.0; + for (auto it = wl.begin(); it != wl.end(); ++it) + sum += *it; + return wl.begin()[k] / sum; + } + } + inline double geometric_pdf(int k, double p) { if (k < 0) - return 0; + return 0.0; else if (k == 0) return p; else return p * std::pow(1 - p, k); } + +#ifdef _GLIBCXX_USE_C99_MATH_TR1 + inline double + negative_binomial_pdf(int k, int n, double p) + { + if (k < 0) + return 0.0; + else + { + double f = std::lgamma(k + (double)n); + double a = std::lgamma(n); + double b = std::lgamma(k + 1.0); + + return std::exp(f - a - b) * std::pow(p, n) * std::pow(1 - p, k); + } + } + + inline double + poisson_pdf(int k, double mu) + { + if (k < 0) + return 0.0; + else + { + double lf = std::lgamma(k + 1.0); + return std::exp(std::log(mu) * k - lf - mu); + } + } +#endif + + inline double + uniform_int_pdf(int k, int a, int b) + { + if (k < 0 || k < a || k > b) + return 0.0; + else + return 1.0 / (b - a + 1.0); + } + } // namespace __gnu_test #endif // #ifndef _GLIBCXX_TESTSUITE_RANDOM_H