Re: more distributions

Re: more distributions

[3dw4rd]

On 03/01/13, Ulrich Drepper<drepper@gmail.com> wrote:
 
I have a few more distributions to be added. The triangle
distribution is the result of combining to uniform distributions and
therefore quite frequently used. The von Mises distribution (the
simple, 2D version) would be the first circular distribution.

The patch depends on the __math_constants patch.

Excellent!

I was looking at a paper: Modeling Data using Directional Distributions by Inderjit S. Dhillon and Suvrit Sra that looks like it would be very similar to your multi-variate normal distribution. These generalize von Mises to higher dimension. Is this your next target?

Thanks for the math constants too.

Ed
 

Re: more distributions

[Ulrich Drepper]

On Fri, Mar 1, 2013 at 5:55 PM,  <3dw4rd@verizon.net> wrote:
> I was looking at a paper: Modeling Data using Directional Distributions by Inderjit S. Dhillon and Suvrit Sra that looks
> like it would be very similar to your multi-variate normal distribution. These generalize von Mises to higher dimension.
> Is this your next target?

von Mises-Fisher?  I might indeed look at it.  The relationship to the
MV Normal distribution is only marginal, though.  MV Normal is not a
circular distribution.

I have a whole bunch of optimization patches which aren't quite done
yet.  I might pick those up first again.

Re: more distributions

[Ulrich Drepper]

On Sat, Mar 2, 2013 at 3:43 PM, Paolo Carlini <paolo.carlini@oracle.com> wrote:
> Yes. Personally, I'm also eager to see your further performance improvements, but I'm afraid will have to wait for 4.9.0.

I checked in the code.  The performance improvements need some
discussions.  I need to do some more experimentation and then will
post some results.

Re: more distributions

[Paolo Carlini]

Hi,

> On Sat, Mar 2, 2013 at 5:21 AM, Paolo Carlini <paolo.carlini@oracle.com> wrote:
>> Exceptionally, I think we can ho ahead with this one too.
> 
> Shall I check in the two patches?

Yes. Personally, I'm also eager to see your further performance improvements, but I'm afraid will have to wait for 4.9.0.

Thanks,
Paolo

Re: more distributions

[Ulrich Drepper]

On Sat, Mar 2, 2013 at 5:21 AM, Paolo Carlini <paolo.carlini@oracle.com> wrote:
> Exceptionally, I think we can ho ahead with this one too.

Shall I check in the two patches?  I added the work-around for
copysign which is the only function used other than log, cos, acos,
sqrt.

Re: more distributions

[Paolo Carlini]

Hi,

On 03/01/2013 07:52 PM, Ulrich Drepper wrote:
> I have a few more distributions to be added.  The triangle
> distribution is the result of combining to uniform distributions and
> therefore quite frequently used.  The von Mises distribution (the
> simple, 2D version) would be the first circular distribution.
Exceptionally, I think we can ho ahead with this one too... *BUT* please 
double check that you are not unconditionally using C99 functions, 
that's a very effective way to break the build on less widespread 
targets. For example I see a std::copysign which needs a replacement 
conditional on _GLIBCXX_USE_C99_MATH_TR1

Paolo.

more distributions

[Ulrich Drepper]

[-- Attachment #1: Type: text/plain, Size: 306 bytes --]

I have a few more distributions to be added.  The triangle
distribution is the result of combining to uniform distributions and
therefore quite frequently used.  The von Mises distribution (the
simple, 2D version) would be the first circular distribution.

The patch depends on the __math_constants patch.

[-- Attachment #2: d-gcc-tri-mises --]
[-- Type: application/octet-stream, Size: 20321 bytes --]

Index: libstdc++-v3/include/ext/random
===================================================================
--- libstdc++-v3/include/ext/random	(revision 196387)
+++ libstdc++-v3/include/ext/random	(working copy)
@@ -37,6 +37,7 @@
 
 #include <random>
 #include <array>
+#include <ext/cmath>
 #ifdef __SSE2__
 # include <x86intrin.h>
 #endif
@@ -2322,6 +2323,513 @@
 	       const hoyt_distribution<_RealType>& __d2)
     { return !(__d1 == __d2); }
 
+
+  /**
+   * @brief A triangular distribution for random numbers.
+   *
+   * The formula for the triangular probability density function is
+   * @f[
+   *                  / 0                          for x < a
+   *     p(x|a,b,c) = | \frac{2(x-a)}{(c-a)(b-a)}  for a <= x <= b
+   *                  | \frac{2(c-x)}{(c-a)(c-b)}  for b < x <= c
+   *                  \ 0                          for c < x
+   * @f]
+   *
+   * <table border=1 cellpadding=10 cellspacing=0>
+   * <caption align=top>Distribution Statistics</caption>
+   * <tr><td>Mean</td><td>@f$ \frac{a+b+c}{2} @f$</td></tr>
+   * <tr><td>Variance</td><td>@f$ \frac{a^2+b^2+c^2-ab-ac-bc}
+   *                                   {18}@f$</td></tr>
+   * <tr><td>Range</td><td>@f$[a, c]@f$</td></tr>
+   * </table>
+   */
+  template<typename _RealType = double>
+    class triangular_distribution
+    {
+      static_assert(std::is_floating_point<_RealType>::value,
+		    "template argument not a floating point type");
+
+    public:
+      /** The type of the range of the distribution. */
+      typedef _RealType result_type;
+      /** Parameter type. */
+      struct param_type
+      {
+	friend class triangular_distribution<_RealType>;
+
+	explicit
+	param_type(_RealType __a = _RealType(0),
+		   _RealType __b = _RealType(0.5),
+		   _RealType __c = _RealType(1))
+	: _M_a(__a), _M_b(__b), _M_c(__c)
+	{
+	  _GLIBCXX_DEBUG_ASSERT(_M_a <= _M_b);
+	  _GLIBCXX_DEBUG_ASSERT(_M_b <= _M_c);
+	  _GLIBCXX_DEBUG_ASSERT(_M_a < _M_c);
+
+	  _M_r_ab = (_M_b - _M_a) / (_M_c - _M_a);
+	  _M_f_ab_ac = (_M_b - _M_a) * (_M_c - _M_a);
+	  _M_f_bc_ac = (_M_c - _M_b) * (_M_c - _M_a);
+	}
+
+	_RealType
+	a() const
+	{ return _M_a; }
+
+	_RealType
+	b() const
+	{ return _M_b; }
+
+	_RealType
+	c() const
+	{ return _M_c; }
+
+	friend bool
+	operator==(const param_type& __p1, const param_type& __p2)
+	{ return (__p1._M_a == __p2._M_a && __p1._M_b == __p2._M_b
+		  && __p1._M_c == __p2._M_c); }
+
+      private:
+
+	_RealType _M_a;
+	_RealType _M_b;
+	_RealType _M_c;
+	_RealType _M_r_ab;
+	_RealType _M_f_ab_ac;
+	_RealType _M_f_bc_ac;
+      };
+
+      /**
+       * @brief Constructs a triangle distribution with parameters
+       * @f$ a @f$, @f$ b @f$ and @f$ c @f$.
+       */
+      explicit
+      triangular_distribution(_RealType __a = _RealType(0),
+			      _RealType __b = _RealType(0.5),
+			      _RealType __c = _RealType(1))
+	: _M_param(__a, __b, __c)
+      { }
+
+      explicit
+      triangular_distribution(const param_type& __p)
+      : _M_param(__p)
+      { }
+
+      /**
+       * @brief Resets the distribution state.
+       */
+      void
+      reset()
+      { }
+
+      /**
+       * @brief Returns the @f$ a @f$ of the distribution.
+       */
+      _RealType
+      a() const
+      { return _M_param.a(); }
+
+      /**
+       * @brief Returns the @f$ b @f$ of the distribution.
+       */
+      _RealType
+      b() const
+      { return _M_param.b(); }
+
+      /**
+       * @brief Returns the @f$ c @f$ of the distribution.
+       */
+      _RealType
+      c() const
+      { return _M_param.c(); }
+
+      /**
+       * @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 _M_param._M_a; }
+
+      /**
+       * @brief Returns the least upper bound value of the distribution.
+       */
+      result_type
+      max() const
+      { return _M_param._M_c; }
+
+      /**
+       * @brief Generating functions.
+       */
+      template<typename _UniformRandomNumberGenerator>
+	result_type
+	operator()(_UniformRandomNumberGenerator& __urng)
+	{ return this->operator()(__urng, _M_param); }
+
+      template<typename _UniformRandomNumberGenerator>
+	result_type
+	operator()(_UniformRandomNumberGenerator& __urng,
+		   const param_type& __p)
+	{
+	  std::__detail::_Adaptor<_UniformRandomNumberGenerator, result_type>
+	    __aurng(__urng);
+	  result_type __rnd = __aurng();
+	  if (__rnd <= __p._M_r_ab)
+	    return __p.a() + std::sqrt(__rnd * __p._M_f_ab_ac);
+	  else
+	    return __p.c() - std::sqrt((result_type(1) - __rnd)
+				       * __p._M_f_bc_ac);
+	}
+
+      template<typename _ForwardIterator,
+	       typename _UniformRandomNumberGenerator>
+	void
+	__generate(_ForwardIterator __f, _ForwardIterator __t,
+		   _UniformRandomNumberGenerator& __urng)
+	{ this->__generate(__f, __t, __urng, _M_param); }
+
+      template<typename _ForwardIterator,
+	       typename _UniformRandomNumberGenerator>
+	void
+	__generate(_ForwardIterator __f, _ForwardIterator __t,
+		   _UniformRandomNumberGenerator& __urng,
+		   const param_type& __p)
+	{ this->__generate_impl(__f, __t, __urng, __p); }
+
+      template<typename _UniformRandomNumberGenerator>
+	void
+	__generate(result_type* __f, result_type* __t,
+		   _UniformRandomNumberGenerator& __urng,
+		   const param_type& __p)
+	{ this->__generate_impl(__f, __t, __urng, __p); }
+
+      /**
+       * @brief Return true if two triangle distributions have the same
+       *        parameters and the sequences that would be generated
+       *        are equal.
+       */
+      friend bool
+      operator==(const triangular_distribution& __d1,
+		 const triangular_distribution& __d2)
+      { return __d1._M_param == __d2._M_param; }
+
+      /**
+       * @brief Inserts a %triangular_distribution random number distribution
+       * @p __x into the output stream @p __os.
+       *
+       * @param __os An output stream.
+       * @param __x  A %triangular_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 __gnu_cxx::triangular_distribution<_RealType1>& __x);
+
+      /**
+       * @brief Extracts a %triangular_distribution random number distribution
+       * @p __x from the input stream @p __is.
+       *
+       * @param __is An input stream.
+       * @param __x  A %triangular_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,
+		   __gnu_cxx::triangular_distribution<_RealType1>& __x);
+
+    private:
+      template<typename _ForwardIterator,
+	       typename _UniformRandomNumberGenerator>
+	void
+	__generate_impl(_ForwardIterator __f, _ForwardIterator __t,
+			_UniformRandomNumberGenerator& __urng,
+			const param_type& __p);
+
+      param_type _M_param;
+    };
+
+  /**
+   * @brief Return true if two triangle distributions are different.
+   */
+  template<typename _RealType>
+    inline bool
+    operator!=(const __gnu_cxx::triangular_distribution<_RealType>& __d1,
+	       const __gnu_cxx::triangular_distribution<_RealType>& __d2)
+   { return !(__d1 == __d2); }
+
+
+  /**
+   * @brief A von Mises distribution for random numbers.
+   *
+   * The formula for the von Mises probability density function is
+   * @f[
+   *     p(x|\mu,\kappa) = \frac{e^{\kappa \cos(x-\mu)}}
+   *                            {2\pi I_0(\kappa)}
+   * @f]
+   *
+   * <table border=1 cellpadding=10 cellspacing=0>
+   * <caption align=top>Distribution Statistics</caption>
+   * <tr><td>Mean</td><td>@f$ \mu @f$</td></tr>
+   * <tr><td>Variance</td><td>@f$ 1-I_1(\kappa)/I_0(\kappa) @f$</td></tr>
+   * <tr><td>Range</td><td>@f$[-\pi, \pi]@f$</td></tr>
+   * </table>
+   */
+  template<typename _RealType = double>
+    class vonmises_distribution
+    {
+      static_assert(std::is_floating_point<_RealType>::value,
+		    "template argument not a floating point type");
+
+    public:
+      /** The type of the range of the distribution. */
+      typedef _RealType result_type;
+      /** Parameter type. */
+      struct param_type
+      {
+	friend class vonmises_distribution<_RealType>;
+
+	explicit
+	param_type(_RealType __mu = _RealType(0),
+		   _RealType __kappa = _RealType(1))
+	: _M_mu(__mu), _M_kappa(__kappa)
+	{
+	  _GLIBCXX_DEBUG_ASSERT(_M_kappa >= _RealType(0));
+	}
+
+	_RealType
+	mu() const
+	{ return _M_mu; }
+
+	_RealType
+	kappa() const
+	{ return _M_kappa; }
+
+	friend bool
+	operator==(const param_type& __p1, const param_type& __p2)
+	{ return __p1._M_kappa == __p2._M_kappa; }
+
+      private:
+
+	_RealType _M_mu;
+	_RealType _M_kappa;
+      };
+
+      /**
+       * @brief Constructs a beta distribution with parameters
+       * @f$\mu@f$ and @f$\kappa@f$.
+       */
+      explicit
+      vonmises_distribution(_RealType __mu = _RealType(0),
+			    _RealType __kappa = _RealType(1))
+	: _M_param(__mu, __kappa)
+      { }
+
+      explicit
+      vonmises_distribution(const param_type& __p)
+      : _M_param(__p)
+      { }
+
+      /**
+       * @brief Resets the distribution state.
+       */
+      void
+      reset()
+      { }
+
+      /**
+       * @brief Returns the @f$ \mu @f$ of the distribution.
+       */
+      _RealType
+      mu() const
+      { return _M_param.mu(); }
+
+      /**
+       * @brief Returns the @f$ \kappa @f$ of the distribution.
+       */
+      _RealType
+      kappa() const
+      { return _M_param.kappa(); }
+
+      /**
+       * @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 this->mu() - __gnu_cxx::__math_constants<_RealType>::__pi;
+      }
+
+      /**
+       * @brief Returns the least upper bound value of the distribution.
+       */
+      result_type
+      max() const
+      {
+	return this->mu() + __gnu_cxx::__math_constants<_RealType>::__pi;
+      }
+
+      /**
+       * @brief Generating functions.
+       */
+      template<typename _UniformRandomNumberGenerator>
+	result_type
+	operator()(_UniformRandomNumberGenerator& __urng)
+	{ return this->operator()(__urng, _M_param); }
+
+      template<typename _UniformRandomNumberGenerator>
+	result_type
+	operator()(_UniformRandomNumberGenerator& __urng,
+		   const param_type& __p)
+	{
+	  const _RealType __pi = __gnu_cxx::__math_constants<_RealType>::__pi;
+	  std::__detail::_Adaptor<_UniformRandomNumberGenerator, result_type>
+	    __aurng(__urng);
+	  result_type __tau = (std::sqrt(result_type(4) * this->kappa()
+					 * this->kappa() + result_type(1))
+			       + result_type(1));
+	  result_type __rho = ((__tau - std::sqrt(result_type(2) * __tau))
+			       / (result_type(2) * this->kappa()));
+	  result_type __r = ((result_type(1) + __rho * __rho)
+			     / (result_type(2) * __rho));
+
+	  result_type __f;
+	  while (1)
+	    {
+	      result_type __rnd = std::cos(__pi * __aurng());
+	      __f = (result_type(1) + __r * __rnd) / (__r + __rnd);
+	      result_type __c = this->kappa() * (__r - __f);
+
+	      result_type __rnd2 = __aurng();
+	      if (__c * (result_type(2) - __c) > __rnd2)
+		break;
+	      if (std::log(__c / __rnd2) >= __c - result_type(1))
+		break;
+	    }
+
+	  result_type __res = std::copysign(std::acos(__f),
+					    __aurng() - result_type(0.5));
+	  if (__res > __pi)
+	    __res -= result_type(2) * __pi;
+	  else if (__res < -__pi)
+	    __res += result_type(2) * __pi;
+	  return __res;
+	}
+
+      template<typename _ForwardIterator,
+	       typename _UniformRandomNumberGenerator>
+	void
+	__generate(_ForwardIterator __f, _ForwardIterator __t,
+		   _UniformRandomNumberGenerator& __urng)
+	{ this->__generate(__f, __t, __urng, _M_param); }
+
+      template<typename _ForwardIterator,
+	       typename _UniformRandomNumberGenerator>
+	void
+	__generate(_ForwardIterator __f, _ForwardIterator __t,
+		   _UniformRandomNumberGenerator& __urng,
+		   const param_type& __p)
+	{ this->__generate_impl(__f, __t, __urng, __p); }
+
+      template<typename _UniformRandomNumberGenerator>
+	void
+	__generate(result_type* __f, result_type* __t,
+		   _UniformRandomNumberGenerator& __urng,
+		   const param_type& __p)
+	{ this->__generate_impl(__f, __t, __urng, __p); }
+
+      /**
+       * @brief Return true if two von Mises distributions have the same
+       *        parameters and the sequences that would be generated
+       *        are equal.
+       */
+      friend bool
+      operator==(const vonmises_distribution& __d1,
+		 const vonmises_distribution& __d2)
+      { return __d1._M_param == __d2._M_param; }
+
+      /**
+       * @brief Inserts a %vonmises_distribution random number distribution
+       * @p __x into the output stream @p __os.
+       *
+       * @param __os An output stream.
+       * @param __x  A %vonmises_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 __gnu_cxx::vonmises_distribution<_RealType1>& __x);
+
+      /**
+       * @brief Extracts a %vonmises_distribution random number distribution
+       * @p __x from the input stream @p __is.
+       *
+       * @param __is An input stream.
+       * @param __x  A %vonmises_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,
+		   __gnu_cxx::vonmises_distribution<_RealType1>& __x);
+
+    private:
+      template<typename _ForwardIterator,
+	       typename _UniformRandomNumberGenerator>
+	void
+	__generate_impl(_ForwardIterator __f, _ForwardIterator __t,
+			_UniformRandomNumberGenerator& __urng,
+			const param_type& __p);
+
+      param_type _M_param;
+    };
+
+  /**
+   * @brief Return true if two von Mises distributions are different.
+   */
+  template<typename _RealType>
+    inline bool
+    operator!=(const __gnu_cxx::vonmises_distribution<_RealType>& __d1,
+	       const __gnu_cxx::vonmises_distribution<_RealType>& __d2)
+   { return !(__d1 == __d2); }
+
 _GLIBCXX_END_NAMESPACE_VERSION
 } // namespace __gnu_cxx
 
Index: libstdc++-v3/include/ext/random.tcc
===================================================================
--- libstdc++-v3/include/ext/random.tcc	(revision 196387)
+++ libstdc++-v3/include/ext/random.tcc	(working copy)
@@ -1187,6 +1187,126 @@
       return __is;
     }
 
+
+  template<typename _RealType>
+    template<typename _OutputIterator,
+	     typename _UniformRandomNumberGenerator>
+      void
+      triangular_distribution<_RealType>::
+      __generate_impl(_OutputIterator __f, _OutputIterator __t,
+		      _UniformRandomNumberGenerator& __urng,
+		      const param_type& __param)
+      {
+	__glibcxx_function_requires(_OutputIteratorConcept<_OutputIterator>)
+
+	while (__f != __t)
+	  *__f++ = this->operator()(__urng, __param);
+      }
+
+  template<typename _RealType, typename _CharT, typename _Traits>
+    std::basic_ostream<_CharT, _Traits>&
+    operator<<(std::basic_ostream<_CharT, _Traits>& __os,
+	       const __gnu_cxx::triangular_distribution<_RealType>& __x)
+    {
+      typedef std::basic_ostream<_CharT, _Traits>  __ostream_type;
+      typedef typename __ostream_type::ios_base    __ios_base;
+
+      const typename __ios_base::fmtflags __flags = __os.flags();
+      const _CharT __fill = __os.fill();
+      const std::streamsize __precision = __os.precision();
+      const _CharT __space = __os.widen(' ');
+      __os.flags(__ios_base::scientific | __ios_base::left);
+      __os.fill(__space);
+      __os.precision(std::numeric_limits<_RealType>::max_digits10);
+
+      __os << __x.a() << __space << __x.b() << __space << __x.c();
+
+      __os.flags(__flags);
+      __os.fill(__fill);
+      __os.precision(__precision);
+      return __os;
+    }
+
+  template<typename _RealType, typename _CharT, typename _Traits>
+    std::basic_istream<_CharT, _Traits>&
+    operator>>(std::basic_istream<_CharT, _Traits>& __is,
+	       __gnu_cxx::triangular_distribution<_RealType>& __x)
+    {
+      typedef std::basic_istream<_CharT, _Traits>  __istream_type;
+      typedef typename __istream_type::ios_base    __ios_base;
+
+      const typename __ios_base::fmtflags __flags = __is.flags();
+      __is.flags(__ios_base::dec | __ios_base::skipws);
+
+      _RealType __a, __b, __c;
+      __is >> __a >> __b >> __c;
+      __x.param(typename __gnu_cxx::triangular_distribution<_RealType>::
+		param_type(__a, __b, __c));
+
+      __is.flags(__flags);
+      return __is;
+    }
+
+
+  template<typename _RealType>
+    template<typename _OutputIterator,
+	     typename _UniformRandomNumberGenerator>
+      void
+      vonmises_distribution<_RealType>::
+      __generate_impl(_OutputIterator __f, _OutputIterator __t,
+		      _UniformRandomNumberGenerator& __urng,
+		      const param_type& __param)
+      {
+	__glibcxx_function_requires(_OutputIteratorConcept<_OutputIterator>)
+
+	while (__f != __t)
+	  *__f++ = this->operator()(__urng, __param);
+      }
+
+  template<typename _RealType, typename _CharT, typename _Traits>
+    std::basic_ostream<_CharT, _Traits>&
+    operator<<(std::basic_ostream<_CharT, _Traits>& __os,
+	       const __gnu_cxx::vonmises_distribution<_RealType>& __x)
+    {
+      typedef std::basic_ostream<_CharT, _Traits>  __ostream_type;
+      typedef typename __ostream_type::ios_base    __ios_base;
+
+      const typename __ios_base::fmtflags __flags = __os.flags();
+      const _CharT __fill = __os.fill();
+      const std::streamsize __precision = __os.precision();
+      const _CharT __space = __os.widen(' ');
+      __os.flags(__ios_base::scientific | __ios_base::left);
+      __os.fill(__space);
+      __os.precision(std::numeric_limits<_RealType>::max_digits10);
+
+      __os << __x.mu() << __space << __x.kappa();
+
+      __os.flags(__flags);
+      __os.fill(__fill);
+      __os.precision(__precision);
+      return __os;
+    }
+
+  template<typename _RealType, typename _CharT, typename _Traits>
+    std::basic_istream<_CharT, _Traits>&
+    operator>>(std::basic_istream<_CharT, _Traits>& __is,
+	       __gnu_cxx::vonmises_distribution<_RealType>& __x)
+    {
+      typedef std::basic_istream<_CharT, _Traits>  __istream_type;
+      typedef typename __istream_type::ios_base    __ios_base;
+
+      const typename __ios_base::fmtflags __flags = __is.flags();
+      __is.flags(__ios_base::dec | __ios_base::skipws);
+
+      _RealType __mu, __kappa;
+      __is >> __mu >> __kappa;
+      __x.param(typename __gnu_cxx::vonmises_distribution<_RealType>::
+		param_type(__mu, __kappa));
+
+      __is.flags(__flags);
+      return __is;
+    }
+
 _GLIBCXX_END_NAMESPACE_VERSION
 } // namespace