Marc Glisse <marc.gli...@inria.fr> writes:

> On Sat, 9 Aug 2014, Ulrich Drepper wrote:
> Yes, you still need the normalization step (divide by the norm).

I guess we can do this.

How about the patch below?  Instead of specializing the entire class for
_Dimen==2 I've added a class at the implementation level.

I've also improved existing tests and add some new ones.


2014-08-09  Ulrich Drepper  <drep...@gmail.com>

        * include/ext/random.tcc (uniform_on_sphere_helper): Define.
        (uniform_on_sphere_distribution::operator()): Use the new helper
        class for the implementation.

        * testsuite/ext/random/uniform_on_sphere_distribution/operators/
        equal.cc: Remove bogus part of comment.
        * testsuite/ext/random/uniform_on_sphere_distribution/operators/
        inequal.cc: Likewise.
        * testsuite/ext/random/uniform_on_sphere_distribution/operators/
        serialize.cc: Add check to verify result of serialzation and
        deserialization.
        * testsuite/ext/random/uniform_on_sphere_distribution/operators/
        generate.cc: New file.


diff --git a/libstdc++-v3/include/ext/random.tcc 
b/libstdc++-v3/include/ext/random.tcc
index 05361d8..d536ecb 100644
--- a/libstdc++-v3/include/ext/random.tcc
+++ b/libstdc++-v3/include/ext/random.tcc
@@ -1540,6 +1540,83 @@ _GLIBCXX_BEGIN_NAMESPACE_VERSION
     }
 
 
+  namespace {
+
+    // Helper class for the uniform_on_sphere_distribution generation
+    // function.
+    template<std::size_t _Dimen, typename _RealType>
+      class uniform_on_sphere_helper
+      {
+       typedef typename uniform_on_sphere_distribution<_Dimen, 
_RealType>::result_type result_type;
+
+      public:
+       template<typename _NormalDistribution, typename 
_UniformRandomNumberGenerator>
+       result_type operator()(_NormalDistribution& __nd,
+                              _UniformRandomNumberGenerator& __urng)
+        {
+         result_type __ret;
+         typename result_type::value_type __norm;
+
+         do
+           {
+             auto __sum = _RealType(0);
+
+             std::generate(__ret.begin(), __ret.end(),
+                           [&__nd, &__urng, &__sum](){
+                             _RealType __t = __nd(__urng);
+                             __sum += __t * __t;
+                             return __t; });
+             __norm = std::sqrt(__sum);
+           }
+         while (__norm == _RealType(0) || ! std::isfinite(__norm));
+
+         std::transform(__ret.begin(), __ret.end(), __ret.begin(),
+                        [__norm](_RealType __val){ return __val / __norm; });
+
+         return __ret;
+        }
+      };
+
+
+    template<typename _RealType>
+      class uniform_on_sphere_helper<2, _RealType>
+      {
+       typedef typename uniform_on_sphere_distribution<2, _RealType>::
+         result_type result_type;
+
+      public:
+       template<typename _NormalDistribution,
+                typename _UniformRandomNumberGenerator>
+       result_type operator()(_NormalDistribution&,
+                              _UniformRandomNumberGenerator& __urng)
+        {
+         result_type __ret;
+         _RealType __sq;
+         std::__detail::_Adaptor<_UniformRandomNumberGenerator,
+                                 _RealType> __aurng(__urng);
+
+         do
+           {
+             __ret[0] = __aurng();
+             __ret[1] = __aurng();
+
+             __sq = __ret[0] * __ret[0] + __ret[1] * __ret[1];
+           }
+         while (__sq == _RealType(0) || __sq > _RealType(1));
+
+         // Yes, we do not just use sqrt(__sq) because hypot() is more
+         // accurate.
+         auto __norm = std::hypot(__ret[0], __ret[1]);
+         __ret[0] /= __norm;
+         __ret[1] /= __norm;
+
+         return __ret;
+        }
+      };
+
+  }
+
+
   template<std::size_t _Dimen, typename _RealType>
     template<typename _UniformRandomNumberGenerator>
       typename uniform_on_sphere_distribution<_Dimen, _RealType>::result_type
@@ -1547,18 +1624,8 @@ _GLIBCXX_BEGIN_NAMESPACE_VERSION
       operator()(_UniformRandomNumberGenerator& __urng,
                 const param_type& __p)
       {
-       result_type __ret;
-       _RealType __sum = _RealType(0);
-
-       std::generate(__ret.begin(), __ret.end(),
-                     [&__urng, &__sum, this](){ _RealType __t = _M_nd(__urng);
-                                                __sum += __t * __t;
-                                                return __t; });
-       auto __norm = std::sqrt(__sum);
-       std::transform(__ret.begin(), __ret.end(), __ret.begin(),
-                      [__norm](_RealType __val){ return __val / __norm; });
-
-       return __ret;
+        uniform_on_sphere_helper<_Dimen, _RealType> __helper;
+        return __helper(_M_nd, __urng);
       }
 
   template<std::size_t _Dimen, typename _RealType>
diff --git 
a/libstdc++-v3/testsuite/ext/random/uniform_on_sphere_distribution/operators/equal.cc
 
b/libstdc++-v3/testsuite/ext/random/uniform_on_sphere_distribution/operators/equal.cc
index 35a024e..f5b8d17 100644
--- 
a/libstdc++-v3/testsuite/ext/random/uniform_on_sphere_distribution/operators/equal.cc
+++ 
b/libstdc++-v3/testsuite/ext/random/uniform_on_sphere_distribution/operators/equal.cc
@@ -20,7 +20,7 @@
 // with this library; see the file COPYING3.  If not see
 // <http://www.gnu.org/licenses/>.
 
-// 26.5.8.4.5 Class template uniform_on_sphere_distribution 
[rand.dist.ext.uniform_on_sphere]
+// Class template uniform_on_sphere_distribution 
[rand.dist.ext.uniform_on_sphere]
 
 #include <ext/random>
 #include <testsuite_hooks.h>
diff --git 
a/libstdc++-v3/testsuite/ext/random/uniform_on_sphere_distribution/operators/inequal.cc
 
b/libstdc++-v3/testsuite/ext/random/uniform_on_sphere_distribution/operators/inequal.cc
index 9f8e8c8..2675652 100644
--- 
a/libstdc++-v3/testsuite/ext/random/uniform_on_sphere_distribution/operators/inequal.cc
+++ 
b/libstdc++-v3/testsuite/ext/random/uniform_on_sphere_distribution/operators/inequal.cc
@@ -20,7 +20,7 @@
 // with this library; see the file COPYING3.  If not see
 // <http://www.gnu.org/licenses/>.
 
-// 26.5.8.4.5 Class template uniform_on_sphere_distribution 
[rand.dist.ext.uniform_on_sphere]
+// Class template uniform_on_sphere_distribution 
[rand.dist.ext.uniform_on_sphere]
 
 #include <ext/random>
 #include <testsuite_hooks.h>
diff --git 
a/libstdc++-v3/testsuite/ext/random/uniform_on_sphere_distribution/operators/serialize.cc
 
b/libstdc++-v3/testsuite/ext/random/uniform_on_sphere_distribution/operators/serialize.cc
index 80264ff..e9a758c 100644
--- 
a/libstdc++-v3/testsuite/ext/random/uniform_on_sphere_distribution/operators/serialize.cc
+++ 
b/libstdc++-v3/testsuite/ext/random/uniform_on_sphere_distribution/operators/serialize.cc
@@ -20,8 +20,8 @@
 // with this library; see the file COPYING3.  If not see
 // <http://www.gnu.org/licenses/>.
 
-// 26.4.8.3.* Class template uniform_on_sphere_distribution 
[rand.dist.ext.uniform_on_sphere]
-// 26.4.2.4 Concept RandomNumberDistribution [rand.concept.dist]
+// Class template uniform_on_sphere_distribution 
[rand.dist.ext.uniform_on_sphere]
+// Concept RandomNumberDistribution [rand.concept.dist]
 
 #include <ext/random>
 #include <sstream>
@@ -40,6 +40,8 @@ test01()
   str << u;
 
   str >> v;
+
+  VERIFY( u == v );
 }
 
 int
--- /dev/null   2014-07-15 08:50:39.432896277 -0400
+++ 
b/libstdc++-v3/testsuite/ext/random/uniform_on_sphere_distribution/operators/generate.cc
    2014-08-09 08:07:29.787244291 -0400
@@ -0,0 +1,60 @@
+// { dg-options "-std=gnu++11" }
+// { dg-require-cstdint "" }
+//
+// 2014-08-09  Ulrich Drepper  <drep...@gmail.com>
+//
+// Copyright (C) 2014 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
+// <http://www.gnu.org/licenses/>.
+
+// Class template uniform_on_sphere_distribution 
[rand.dist.ext.uniform_on_sphere]
+// Concept RandomNumberDistribution [rand.concept.dist]
+
+#include <ext/random>
+#include <sstream>
+#include <testsuite_hooks.h>
+
+void
+test01()
+{
+  bool test [[gnu::unused]] = true;
+  std::minstd_rand0 rng;
+
+  __gnu_cxx::uniform_on_sphere_distribution<3> u3;
+
+  for (size_t n = 0; n < 1000; ++n)
+    {
+      auto r = u3(rng);
+
+      VERIFY (r[0] != 0.0 || r[1] != 0.0 || r[2] != 0.0);
+    }
+
+  __gnu_cxx::uniform_on_sphere_distribution<2> u2;
+
+  for (size_t n = 0; n < 1000; ++n)
+    {
+      auto r = u2(rng);
+
+      VERIFY (r[0] != 0.0 || r[1] != 0.0);
+    }
+}
+
+int
+main()
+{
+  test01();
+  return 0;
+}

Reply via email to