+2006-07-06 Paolo Carlini <pcarlini@suse.de>
+
+ * include/tr1/random (class gamma_distribution<>): Add.
+ * include/tr1/random.tcc (gamma_distribution<>::operator(),
+ operator<<(std::basic_ostream<>&, const gamma_distribution<>&)):
+ Define.
+ * testsuite/tr1/5_numerical_facilities/random/gamma_distribution/
+ requirements/typedefs.cc: New.
+
2006-07-06 Benjamin Kosnik <bkoz@redhat.com>
* testsuite/util/regression/trait/assoc/trait.hpp: Format.
bool _M_saved_available;
};
+
+ /**
+ * @brief A gamma continuous distribution for random numbers.
+ *
+ * The formula for the gamma probability mass function is
+ * @f$ p(x) = \frac{1}{\Gamma(\alpha)} x^{\alpha - 1} e^{-x} } @f$.
+ */
+ template<typename _RealType = double>
+ class gamma_distribution;
+
+ template<typename _RealType, typename _CharT, typename _Traits>
+ std::basic_ostream<_CharT, _Traits>&
+ operator<<(std::basic_ostream<_CharT, _Traits>& __os,
+ const gamma_distribution<_RealType>& __x);
+
+ template<typename _RealType>
+ class gamma_distribution
+ {
+ public:
+ // types
+ typedef _RealType input_type;
+ typedef _RealType result_type;
+
+ public:
+ /**
+ * Constructs a gamma distribution with parameters @f$ \alpha @f$.
+ */
+ explicit
+ gamma_distribution(const result_type& __alpha = result_type(1))
+ : _M_alpha(__alpha)
+ {
+ _GLIBCXX_DEBUG_ASSERT(_M_alpha > 0);
+ }
+
+ /**
+ * Gets the @f$ \alpha @f$ of the distribution.
+ */
+ _RealType
+ alpha() const
+ { return _M_alpha; }
+
+ /**
+ * Resets the distribution.
+ */
+ void
+ reset() { }
+
+ template<class _UniformRandomNumberGenerator>
+ result_type
+ operator()(_UniformRandomNumberGenerator& __urng);
+
+ /**
+ * 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 _RealType1, typename _CharT, typename _Traits>
+ friend std::basic_ostream<_CharT, _Traits>&
+ operator<<(std::basic_ostream<_CharT, _Traits>& __os,
+ const gamma_distribution<_RealType1>& __x);
+
+ /**
+ * 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 _RealType1, typename _CharT, typename _Traits>
+ friend std::basic_istream<_CharT, _Traits>&
+ operator>>(std::basic_istream<_CharT, _Traits>& __is,
+ gamma_distribution<_RealType1>& __x)
+ { return __is >> __x._M_alpha; }
+
+ private:
+ result_type _M_alpha;
+ };
+
/* @} */ // group tr1_random_distributions_continuous
/* @} */ // group tr1_random_distributions
/* @} */ // group tr1_random
return __is;
}
+
+ /**
+ * Cheng's rejection algorithm GB for alpha >= 1 and a modification
+ * of Vaduva's rejection from Weibull algorithm due to Devroye for
+ * alpha < 1.
+ *
+ * References:
+ * Cheng, R. C. "The Generation of Gamma Random Variables with Non-integral
+ * Shape Parameter." Applied Statistics, 26, 71-75, 1977.
+ *
+ * Vaduva, I. "Computer Generation of Gamma Gandom Variables by Rejection
+ * and Composition Procedures." Math. Operationsforschung and Statistik,
+ * Series in Statistics, 8, 545-576, 1977.
+ *
+ * Devroye, L. "Non-Uniform Random Variates Generation." Springer-Verlag,
+ * New York, 1986, Sect. 3.4.
+ */
+ template<typename _RealType>
+ template<class _UniformRandomNumberGenerator>
+ typename gamma_distribution<_RealType>::result_type
+ gamma_distribution<_RealType>::
+ operator()(_UniformRandomNumberGenerator& __urng)
+ {
+ result_type __x;
+
+ if (_M_alpha >= 1)
+ {
+ // alpha - log(4)
+ const result_type __b = _M_alpha
+ - result_type(1.3862943611198906188344642429163531L);
+ const result_type __c = _M_alpha + std::sqrt(2 * _M_alpha - 1);
+
+ // 1 + log(9 / 2)
+ const result_type __k = 2.5040773967762740733732583523868748L;
+
+ result_type __z, __r;
+ do
+ {
+ const result_type __u = __urng();
+ const result_type __v = __urng();
+
+ const result_type __y = _M_alpha * std::log(__v / (1 - __v));
+ __x = _M_alpha * std::exp(__v);
+
+ __z = __u * __v * __v;
+ __r = __b + __c * __y - __x;
+ }
+ while (__r < result_type(4.5) * __z - __k
+ && __r < std::log(__z));
+ }
+ else
+ {
+ const result_type __c = 1 / _M_alpha;
+ const result_type __d =
+ std::pow(_M_alpha, _M_alpha / (1 - _M_alpha)) * (1 - _M_alpha);
+
+ result_type __z, __e;
+ do
+ {
+ __z = -std::log(__urng());
+ __e = -std::log(__urng());
+
+ __x = std::pow(__z, __c);
+ }
+ while (__z + __e > __d + __x);
+ }
+
+ return __x;
+ }
+
+ template<typename _RealType, typename _CharT, typename _Traits>
+ std::basic_ostream<_CharT, _Traits>&
+ operator<<(std::basic_ostream<_CharT, _Traits>& __os,
+ const gamma_distribution<_RealType>& __x)
+ {
+ const std::ios_base::fmtflags __flags = __os.flags();
+ const _CharT __fill = __os.fill();
+ const std::streamsize __precision = __os.precision();
+ __os.flags(std::ios_base::scientific | std::ios_base::left);
+ __os.fill(__os.widen(' '));
+ __os.precision(_Private::_Max_digits10<_RealType>::__value);
+
+ __os << __x.alpha();
+
+ __os.flags(__flags);
+ __os.fill(__fill);
+ __os.precision(__precision);
+ return __os;
+ }
+
_GLIBCXX_END_NAMESPACE
}
--- /dev/null
+// { dg-do compile }
+//
+// 2006-07-06 Paolo Carlini <pcarlini@suse.de>
+//
+// Copyright (C) 2006 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 2, 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 COPYING. If not, write to the Free
+// Software Foundation, 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301,
+// USA.
+
+// 5.1.7.9 Class template gamma_distribution [tr.rand.dist.gamma]
+// 5.1.1 [7] Table 17
+
+#include <tr1/random>
+
+void
+test01()
+{
+ using namespace std::tr1;
+
+ typedef gamma_distribution<double> test_type;
+
+ typedef test_type::input_type input_type;
+ typedef test_type::result_type result_type;
+}