This is the mail archive of the
gcc-patches@gcc.gnu.org
mailing list for the GCC project.
[v3] Some more work on <random>
- From: Paolo Carlini <paolo dot carlini at oracle dot com>
- To: "gcc-patches at gcc dot gnu dot org" <gcc-patches at gcc dot gnu dot org>
- Cc: libstdc++ at gcc dot gnu dot org
- Date: Thu, 11 Jun 2009 20:34:13 +0200
- Subject: [v3] Some more work on <random>
Hi,
tested x86_64-linux multilib, committed.
Paolo.
///////////////////////
2009-06-11 Paolo Carlini <paolo.carlini@oracle.com>
* include/bits/random.tcc
(negative_binomial_distribution<>::operator()
(_UniformRandomNumberGenerator&, const param_type&): Tweak to use a
class member gamma_distribution.
(negative_binomial_distribution<>::operator()
(_UniformRandomNumberGenerator&)): Implement out of line here.
(operator<<(basic_ostream<>&, negative_binomial_distribution<>),
operator>>(basic_ostream<>&, negative_binomial_distribution<>): Adjust.
(student_t_distribution<>::operator()
(_UniformRandomNumberGenerator&, const param_type&): Move inline,
simplify.
(operator<<(basic_ostream<>&, student_t_distribution<>),
operator>>(basic_ostream<>&, student_t_distribution<>): Adjust.
(chi_squared_distribution<>::operator()
(_UniformRandomNumberGenerator&, const param_type&): Move inline,
tweak to use a class member gamma_distribution.
(operator<<(basic_ostream<>&, chi_squared_distribution<>),
operator>>(basic_ostream<>&, chi_squared_distribution<>): Adjust.
(fisher_f_distribution<>::operator() (_UniformRandomNumberGenerator&,
const param_type&): Move inline, tweak to use class member
gamma_distributions.
(operator<<(basic_ostream<>&, fisher_f_distribution<>),
operator>>(basic_ostream<>&, fisher_f_distribution<>): Adjust.
* include/bits/random.h: Adjust, minor tweaks.
Index: include/bits/random.tcc
===================================================================
--- include/bits/random.tcc (revision 148392)
+++ include/bits/random.tcc (working copy)
@@ -854,18 +854,34 @@
return __is;
}
+
template<typename _IntType>
template<typename _UniformRandomNumberGenerator>
typename negative_binomial_distribution<_IntType>::result_type
negative_binomial_distribution<_IntType>::
+ operator()(_UniformRandomNumberGenerator& __urng)
+ {
+ const double __y = _M_gd(__urng);
+
+ // XXX Is the constructor too slow?
+ std::poisson_distribution<result_type> __poisson(__y);
+ return __poisson(__urng);
+ }
+
+ template<typename _IntType>
+ template<typename _UniformRandomNumberGenerator>
+ typename negative_binomial_distribution<_IntType>::result_type
+ negative_binomial_distribution<_IntType>::
operator()(_UniformRandomNumberGenerator& __urng,
const param_type& __p)
{
- gamma_distribution<> __gamma(__p.k(), 1.0);
- double __x = __gamma(__urng);
+ typedef typename std::gamma_distribution<result_type>::param_type
+ param_type;
+
+ const double __y =
+ _M_gd(__urng, param_type(__p.k(), __p.p() / (1.0 - __p.p())));
- poisson_distribution<result_type> __poisson(__x * __p.p()
- / (1.0 - __p.p()));
+ std::poisson_distribution<result_type> __poisson(__y);
return __poisson(__urng);
}
@@ -885,7 +901,8 @@
__os.fill(__os.widen(' '));
__os.precision(std::numeric_limits<double>::digits10 + 1);
- __os << __x.k() << __space << __x.p();
+ __os << __x.k() << __space << __x.p()
+ << __space << __x._M_gd;
__os.flags(__flags);
__os.fill(__fill);
@@ -906,7 +923,7 @@
_IntType __k;
double __p;
- __is >> __k >> __p;
+ __is >> __k >> __p >> __x._M_gd;
__x.param(typename negative_binomial_distribution<_IntType>::
param_type(__k, __p));
@@ -1538,17 +1555,6 @@
}
- template<typename _RealType>
- template<typename _UniformRandomNumberGenerator>
- typename chi_squared_distribution<_RealType>::result_type
- chi_squared_distribution<_RealType>::
- operator()(_UniformRandomNumberGenerator& __urng,
- const param_type& __p)
- {
- gamma_distribution<_RealType> __gamma(__p.n() / 2, 1.0);
- return 2 * __gamma(__urng);
- }
-
template<typename _RealType, typename _CharT, typename _Traits>
std::basic_ostream<_CharT, _Traits>&
operator<<(std::basic_ostream<_CharT, _Traits>& __os,
@@ -1565,7 +1571,7 @@
__os.fill(__space);
__os.precision(std::numeric_limits<_RealType>::digits10 + 1);
- __os << __x.n();
+ __os << __x.n() << __space << __x._M_gd;
__os.flags(__flags);
__os.fill(__fill);
@@ -1585,7 +1591,7 @@
__is.flags(__ios_base::dec | __ios_base::skipws);
_RealType __n;
- __is >> __n;
+ __is >> __n >> __x._M_gd;
__x.param(typename chi_squared_distribution<_RealType>::
param_type(__n));
@@ -1657,23 +1663,6 @@
}
- template<typename _RealType>
- template<typename _UniformRandomNumberGenerator>
- typename fisher_f_distribution<_RealType>::result_type
- fisher_f_distribution<_RealType>::
- operator()(_UniformRandomNumberGenerator& __urng,
- const param_type& __p)
- {
- gamma_distribution<_RealType> __gamma;
- _RealType __ym = __gamma(__urng,
- typename gamma_distribution<_RealType>::param_type(__p.m() / 2, 2));
-
- _RealType __yn = __gamma(__urng,
- typename gamma_distribution<_RealType>::param_type(__p.n() / 2, 2));
-
- return (__ym * __p.n()) / (__yn * __p.m());
- }
-
template<typename _RealType, typename _CharT, typename _Traits>
std::basic_ostream<_CharT, _Traits>&
operator<<(std::basic_ostream<_CharT, _Traits>& __os,
@@ -1690,7 +1679,8 @@
__os.fill(__space);
__os.precision(std::numeric_limits<_RealType>::digits10 + 1);
- __os << __x.m() << __space << __x.n();
+ __os << __x.m() << __space << __x.n()
+ << __space << __x._M_gd_x << __space << __x._M_gd_y;
__os.flags(__flags);
__os.fill(__fill);
@@ -1710,7 +1700,7 @@
__is.flags(__ios_base::dec | __ios_base::skipws);
_RealType __m, __n;
- __is >> __m >> __n;
+ __is >> __m >> __n >> __x._M_gd_x >> __x._M_gd_y;
__x.param(typename fisher_f_distribution<_RealType>::
param_type(__m, __n));
@@ -1719,43 +1709,6 @@
}
- template<typename _RealType>
- template<typename _UniformRandomNumberGenerator>
- typename student_t_distribution<_RealType>::result_type
- student_t_distribution<_RealType>::
- operator()(_UniformRandomNumberGenerator& __urng,
- const param_type& __param)
- {
- if (__param.n() <= 2.0)
- {
- _RealType __y1 = _M_nd(__urng);
- chi_squared_distribution<_RealType> __chisq(__param.n());
- _RealType __y2 = __chisq(__urng);
-
- return __y1 / std::sqrt(__y2 / __param.n());
- }
- else
- {
- _RealType __y1, __y2, __z;
- exponential_distribution<_RealType>
- __exponential(1.0 / (__param.n() / 2.0 - 1.0));
-
- do
- {
- __y1 = _M_nd(__urng);
- __y2 = __exponential(__urng);
-
- __z = __y1 * __y1 / (__param.n() - 2.0);
- }
- while (1.0 - __z < 0.0 || std::exp(-__y2 - __z) > (1.0 - __z));
-
- // Note that there is a typo in Knuth's formula, the line below
- // is taken from the original paper of Marsaglia, Mathematics of
- // Computation, 34 (1980), p 234-256
- return __y1 / std::sqrt((1.0 - 2.0 / __param.n()) * (1.0 - __z));
- }
- }
-
template<typename _RealType, typename _CharT, typename _Traits>
std::basic_ostream<_CharT, _Traits>&
operator<<(std::basic_ostream<_CharT, _Traits>& __os,
@@ -1772,7 +1725,7 @@
__os.fill(__space);
__os.precision(std::numeric_limits<_RealType>::digits10 + 1);
- __os << __x.n() << __space << __x._M_nd;
+ __os << __x.n() << __space << __x._M_nd << __space << __x._M_gd;
__os.flags(__flags);
__os.fill(__fill);
@@ -1792,7 +1745,7 @@
__is.flags(__ios_base::dec | __ios_base::skipws);
_RealType __n;
- __is >> __n >> __x._M_nd;
+ __is >> __n >> __x._M_nd >> __x._M_gd;
__x.param(typename student_t_distribution<_RealType>::param_type(__n));
__is.flags(__flags);
Index: include/bits/random.h
===================================================================
--- include/bits/random.h (revision 148392)
+++ include/bits/random.h (working copy)
@@ -2078,7 +2078,7 @@
private:
param_type _M_param;
- normal_distribution<result_type> _M_nd;
+ std::normal_distribution<result_type> _M_nd;
};
/**
@@ -2111,8 +2111,166 @@
operator>>(std::basic_istream<_CharT, _Traits>&,
std::lognormal_distribution<_RealType>&);
+
+ /**
+ * @brief A gamma continuous distribution for random numbers.
+ *
+ * The formula for the gamma probability density function is
+ * @f$ p(x|\alpha,\beta) = \frac{1}{\beta\Gamma(\alpha)}
+ * (x/\beta)^{\alpha - 1} e^{-x/\beta} @f$.
+ */
+ template<typename _RealType = double>
+ class gamma_distribution
+ {
+ public:
+ /** The type of the range of the distribution. */
+ typedef _RealType result_type;
+ /** Parameter type. */
+ struct param_type
+ {
+ typedef gamma_distribution<_RealType> distribution_type;
+ friend class gamma_distribution<_RealType>;
+ explicit
+ param_type(_RealType __alpha_val = _RealType(1),
+ _RealType __beta_val = _RealType(1))
+ : _M_alpha(__alpha_val), _M_beta(__beta_val)
+ {
+ _GLIBCXX_DEBUG_ASSERT(_M_alpha > _RealType(0));
+ _M_initialize();
+ }
+
+ _RealType
+ alpha() const
+ { return _M_alpha; }
+
+ _RealType
+ beta() const
+ { return _M_beta; }
+
+ private:
+ void
+ _M_initialize();
+
+ _RealType _M_alpha;
+ _RealType _M_beta;
+
+ _RealType _M_malpha, _M_a2;
+ };
+
+ public:
+ /**
+ * @brief Constructs a gamma distribution with parameters
+ * @f$ \alpha @f$ and @f$ \beta @f$.
+ */
+ explicit
+ gamma_distribution(_RealType __alpha_val = _RealType(1),
+ _RealType __beta_val = _RealType(1))
+ : _M_param(__alpha_val, __beta_val), _M_nd()
+ { }
+
+ explicit
+ gamma_distribution(const param_type& __p)
+ : _M_param(__p), _M_nd()
+ { }
+
+ /**
+ * @brief Resets the distribution state.
+ */
+ void
+ reset()
+ { _M_nd.reset(); }
+
+ /**
+ * @brief Returns the @f$ \alpha @f$ of the distribution.
+ */
+ _RealType
+ alpha() const
+ { return _M_param.alpha(); }
+
+ /**
+ * @brief Returns the @f$ \beta @f$ of the distribution.
+ */
+ _RealType
+ beta() const
+ { return _M_param.beta(); }
+
+ /**
+ * @brief Returns the parameter set of the distribution.
+ */
+ param_type
+ param() const
+ { return _M_param; }
+
+ /**
+ * @brief Sets the parameter set of the distribution.
+ * @param __param The new parameter set of the distribution.
+ */
+ void
+ param(const param_type& __param)
+ { _M_param = __param; }
+
+ /**
+ * @brief Returns the greatest lower bound value of the distribution.
+ */
+ result_type
+ min() const
+ { return result_type(0); }
+
+ /**
+ * @brief Returns the least upper bound value of the distribution.
+ */
+ result_type
+ max() const
+ { return std::numeric_limits<result_type>::max(); }
+
+ template<typename _UniformRandomNumberGenerator>
+ result_type
+ operator()(_UniformRandomNumberGenerator& __urng)
+ { return this->operator()(__urng, this->param()); }
+
+ template<typename _UniformRandomNumberGenerator>
+ result_type
+ operator()(_UniformRandomNumberGenerator& __urng,
+ const param_type& __p);
+
+ private:
+ param_type _M_param;
+
+ std::normal_distribution<result_type> _M_nd;
+ };
+
/**
+ * @brief Inserts a %gamma_distribution random number distribution
+ * @p __x into the output stream @p __os.
+ *
+ * @param __os An output stream.
+ * @param __x A %gamma_distribution random number distribution.
+ *
+ * @returns The output stream with the state of @p __x inserted or in
+ * an error state.
+ */
+ template<typename _RealType, typename _CharT, typename _Traits>
+ std::basic_ostream<_CharT, _Traits>&
+ operator<<(std::basic_ostream<_CharT, _Traits>&,
+ const std::gamma_distribution<_RealType>&);
+
+ /**
+ * @brief Extracts a %gamma_distribution random number distribution
+ * @p __x from the input stream @p __is.
+ *
+ * @param __is An input stream.
+ * @param __x A %gamma_distribution random number generator engine.
+ *
+ * @returns The input stream with @p __x extracted or in an error state.
+ */
+ template<typename _RealType, typename _CharT, typename _Traits>
+ std::basic_istream<_CharT, _Traits>&
+ operator>>(std::basic_istream<_CharT, _Traits>&,
+ std::gamma_distribution<_RealType>&);
+
+
+ /**
* @brief A chi_squared_distribution random number distribution.
*
* The formula for the normal probability mass function is
@@ -2144,12 +2302,12 @@
explicit
chi_squared_distribution(_RealType __n = _RealType(1))
- : _M_param(__n)
+ : _M_param(__n), _M_gd(__n / 2)
{ }
explicit
chi_squared_distribution(const param_type& __p)
- : _M_param(__p)
+ : _M_param(__p), _M_gd(__p.n() / 2)
{ }
/**
@@ -2157,7 +2315,7 @@
*/
void
reset()
- { }
+ { _M_gd.reset(); }
/**
*
@@ -2198,15 +2356,22 @@
template<typename _UniformRandomNumberGenerator>
result_type
operator()(_UniformRandomNumberGenerator& __urng)
- { return this->operator()(__urng, this->param()); }
+ { return 2 * _M_gd(__urng); }
template<typename _UniformRandomNumberGenerator>
result_type
operator()(_UniformRandomNumberGenerator& __urng,
- const param_type& __p);
+ const param_type& __p)
+ {
+ typedef typename std::gamma_distribution<result_type>::param_type
+ param_type;
+ return 2 * _M_gd(__urng, param_type(__p.n() / 2));
+ }
private:
param_type _M_param;
+
+ std::gamma_distribution<result_type> _M_gd;
};
/**
@@ -2420,12 +2585,12 @@
explicit
fisher_f_distribution(_RealType __m = _RealType(1),
_RealType __n = _RealType(1))
- : _M_param(__m, __n)
+ : _M_param(__m, __n), _M_gd_x(__m / 2), _M_gd_y(__n / 2)
{ }
explicit
fisher_f_distribution(const param_type& __p)
- : _M_param(__p)
+ : _M_param(__p), _M_gd_x(__p.m() / 2), _M_gd_y(__p.n() / 2)
{ }
/**
@@ -2433,7 +2598,10 @@
*/
void
reset()
- { }
+ {
+ _M_gd_x.reset();
+ _M_gd_y.reset();
+ }
/**
*
@@ -2478,15 +2646,23 @@
template<typename _UniformRandomNumberGenerator>
result_type
operator()(_UniformRandomNumberGenerator& __urng)
- { return this->operator()(__urng, this->param()); }
+ { return (_M_gd_x(__urng) * n()) / (_M_gd_y(__urng) * m()); }
template<typename _UniformRandomNumberGenerator>
result_type
operator()(_UniformRandomNumberGenerator& __urng,
- const param_type& __p);
+ const param_type& __p)
+ {
+ typedef typename std::gamma_distribution<result_type>::param_type
+ param_type;
+ return ((_M_gd_x(__urng, param_type(__p.m() / 2)) * n())
+ / (_M_gd_y(__urng, param_type(__p.n() / 2)) * m()));
+ }
private:
param_type _M_param;
+
+ std::gamma_distribution<result_type> _M_gd_x, _M_gd_y;
};
/**
@@ -2553,12 +2729,12 @@
explicit
student_t_distribution(_RealType __n = _RealType(1))
- : _M_param(__n), _M_nd()
+ : _M_param(__n), _M_nd(), _M_gd(__n / 2, 2)
{ }
explicit
student_t_distribution(const param_type& __p)
- : _M_param(__p), _M_nd()
+ : _M_param(__p), _M_nd(), _M_gd(__p.n() / 2, 2)
{ }
/**
@@ -2566,7 +2742,10 @@
*/
void
reset()
- { _M_nd.reset(); }
+ {
+ _M_nd.reset();
+ _M_gd.reset();
+ }
/**
*
@@ -2606,18 +2785,26 @@
template<typename _UniformRandomNumberGenerator>
result_type
- operator()(_UniformRandomNumberGenerator& __urng)
- { return this->operator()(__urng, this->param()); }
+ operator()(_UniformRandomNumberGenerator& __urng)
+ { return _M_nd(__urng) * std::sqrt(n() / _M_gd(__urng)); }
template<typename _UniformRandomNumberGenerator>
result_type
operator()(_UniformRandomNumberGenerator& __urng,
- const param_type& __p);
+ const param_type& __p)
+ {
+ typedef typename std::gamma_distribution<result_type>::param_type
+ param_type;
+
+ const result_type __g = _M_gd(__urng, param_type(__p.n() / 2, 2));
+ return _M_nd(__urng) * std::sqrt(__p.n() / __g);
+ }
private:
param_type _M_param;
- normal_distribution<result_type> _M_nd;
+ std::normal_distribution<result_type> _M_nd;
+ std::gamma_distribution<result_type> _M_gd;
};
/**
@@ -2977,7 +3164,7 @@
param_type _M_param;
// NB: Unused when _GLIBCXX_USE_C99_MATH_TR1 is undefined.
- normal_distribution<double> _M_nd;
+ std::normal_distribution<double> _M_nd;
};
@@ -3166,12 +3353,12 @@
explicit
negative_binomial_distribution(_IntType __k = 1, double __p = 0.5)
- : _M_param(__k, __p)
+ : _M_param(__k, __p), _M_gd(__k, __p / (1.0 - __p))
{ }
explicit
negative_binomial_distribution(const param_type& __p)
- : _M_param(__p)
+ : _M_param(__p), _M_gd(__p.k(), __p.p() / (1.0 - __p.p()))
{ }
/**
@@ -3179,7 +3366,7 @@
*/
void
reset()
- { }
+ { _M_gd.reset(); }
/**
* @brief Return the @f$ k @f$ parameter of the distribution.
@@ -3226,8 +3413,7 @@
template<typename _UniformRandomNumberGenerator>
result_type
- operator()(_UniformRandomNumberGenerator& __urng)
- { return this->operator()(__urng, this->param()); }
+ operator()(_UniformRandomNumberGenerator& __urng);
template<typename _UniformRandomNumberGenerator>
result_type
@@ -3236,6 +3422,8 @@
private:
param_type _M_param;
+
+ std::gamma_distribution<double> _M_gd;
};
/**
@@ -3421,7 +3609,7 @@
param_type _M_param;
// NB: Unused when _GLIBCXX_USE_C99_MATH_TR1 is undefined.
- normal_distribution<double> _M_nd;
+ std::normal_distribution<double> _M_nd;
};
/**
@@ -3575,164 +3763,6 @@
/**
- * @brief A gamma continuous distribution for random numbers.
- *
- * The formula for the gamma probability density function is
- * @f$ p(x|\alpha,\beta) = \frac{1}{\beta\Gamma(\alpha)}
- * (x/\beta)^{\alpha - 1} e^{-x/\beta} @f$.
- */
- template<typename _RealType = double>
- class gamma_distribution
- {
- public:
- /** The type of the range of the distribution. */
- typedef _RealType result_type;
- /** Parameter type. */
- struct param_type
- {
- typedef gamma_distribution<_RealType> distribution_type;
- friend class gamma_distribution<_RealType>;
-
- explicit
- param_type(_RealType __alpha_val = _RealType(1),
- _RealType __beta_val = _RealType(1))
- : _M_alpha(__alpha_val), _M_beta(__beta_val)
- {
- _GLIBCXX_DEBUG_ASSERT(_M_alpha > _RealType(0));
- _M_initialize();
- }
-
- _RealType
- alpha() const
- { return _M_alpha; }
-
- _RealType
- beta() const
- { return _M_beta; }
-
- private:
- void
- _M_initialize();
-
- _RealType _M_alpha;
- _RealType _M_beta;
-
- _RealType _M_malpha, _M_a2;
- };
-
- public:
- /**
- * @brief Constructs a gamma distribution with parameters
- * @f$ \alpha @f$ and @f$ \beta @f$.
- */
- explicit
- gamma_distribution(_RealType __alpha_val = _RealType(1),
- _RealType __beta_val = _RealType(1))
- : _M_param(__alpha_val, __beta_val), _M_nd()
- { }
-
- explicit
- gamma_distribution(const param_type& __p)
- : _M_param(__p), _M_nd()
- { }
-
- /**
- * @brief Resets the distribution state.
- */
- void
- reset()
- { _M_nd.reset(); }
-
- /**
- * @brief Returns the @f$ \alpha @f$ of the distribution.
- */
- _RealType
- alpha() const
- { return _M_param.alpha(); }
-
- /**
- * @brief Returns the @f$ \beta @f$ of the distribution.
- */
- _RealType
- beta() const
- { return _M_param.beta(); }
-
- /**
- * @brief Returns the parameter set of the distribution.
- */
- param_type
- param() const
- { return _M_param; }
-
- /**
- * @brief Sets the parameter set of the distribution.
- * @param __param The new parameter set of the distribution.
- */
- void
- param(const param_type& __param)
- { _M_param = __param; }
-
- /**
- * @brief Returns the greatest lower bound value of the distribution.
- */
- result_type
- min() const
- { return result_type(0); }
-
- /**
- * @brief Returns the least upper bound value of the distribution.
- */
- result_type
- max() const
- { return std::numeric_limits<result_type>::max(); }
-
- template<typename _UniformRandomNumberGenerator>
- result_type
- operator()(_UniformRandomNumberGenerator& __urng)
- { return this->operator()(__urng, this->param()); }
-
- template<typename _UniformRandomNumberGenerator>
- result_type
- operator()(_UniformRandomNumberGenerator& __urng,
- const param_type& __p);
-
- private:
- param_type _M_param;
-
- normal_distribution<result_type> _M_nd;
- };
-
- /**
- * @brief Inserts a %gamma_distribution random number distribution
- * @p __x into the output stream @p __os.
- *
- * @param __os An output stream.
- * @param __x A %gamma_distribution random number distribution.
- *
- * @returns The output stream with the state of @p __x inserted or in
- * an error state.
- */
- template<typename _RealType, typename _CharT, typename _Traits>
- std::basic_ostream<_CharT, _Traits>&
- operator<<(std::basic_ostream<_CharT, _Traits>&,
- const std::gamma_distribution<_RealType>&);
-
- /**
- * @brief Extracts a %gamma_distribution random number distribution
- * @p __x from the input stream @p __is.
- *
- * @param __is An input stream.
- * @param __x A %gamma_distribution random number generator engine.
- *
- * @returns The input stream with @p __x extracted or in an error state.
- */
- template<typename _RealType, typename _CharT, typename _Traits>
- std::basic_istream<_CharT, _Traits>&
- operator>>(std::basic_istream<_CharT, _Traits>&,
- std::gamma_distribution<_RealType>&);
-
-
- /**
* @brief A weibull_distribution random number distribution.
*
* The formula for the normal probability density function is