Re: [PATCH] libstdc++: add uniform on sphere distribution
On 08/13/2014 01:32 PM, Marc Glisse wrote: On Wed, 13 Aug 2014, Paolo Carlini wrote: ... fixed with the below. Could you please use __sq instead of recomputing it? Patch preapproved, thanks! Paolo.
Re: [PATCH] libstdc++: add uniform on sphere distribution
On Wed, 13 Aug 2014, Paolo Carlini wrote: ... fixed with the below. Could you please use __sq instead of recomputing it? -- Marc Glisse
Re: [PATCH] libstdc++: add uniform on sphere distribution
... fixed with the below.
Thanks,
Paolo.
//
2014-08-13 Paolo Carlini
PR libstdc++/62118
* include/ext/random.tcc (uniform_on_sphere_helper<2, _RealType>::
operator()): Use std::hypot only when _GLIBCXX_USE_C99_MATH_TR1.
Index: include/ext/random.tcc
===
--- include/ext/random.tcc (revision 213898)
+++ include/ext/random.tcc (working copy)
@@ -1547,10 +1547,12 @@
template
class uniform_on_sphere_helper
{
- typedef typename uniform_on_sphere_distribution<_Dimen,
_RealType>::result_type result_type;
+ typedef typename uniform_on_sphere_distribution<_Dimen, _RealType>::
+ result_type result_type;
public:
- template
+ template
result_type operator()(_NormalDistribution& __nd,
_UniformRandomNumberGenerator& __urng)
{
@@ -1604,9 +1606,13 @@
}
while (__sq == _RealType(0) || __sq > _RealType(1));
+#if _GLIBCXX_USE_C99_MATH_TR1
// Yes, we do not just use sqrt(__sq) because hypot() is more
// accurate.
auto __norm = std::hypot(__ret[0], __ret[1]);
+#else
+ auto __norm = std::sqrt(__ret[0] * __ret[0] + __ret[1] * __ret[1]);
+#endif
__ret[0] /= __norm;
__ret[1] /= __norm;
Re: [PATCH] libstdc++: add uniform on sphere distribution
Ulrich-- On 08/13/2014 09:36 AM, Bin.Cheng wrote: On Sun, Aug 10, 2014 at 2:00 AM, Ulrich Drepper wrote: On Sat, Aug 9, 2014 at 1:40 PM, Marc Glisse wrote: __x = result_type(2.0) * __aurng() - 1.0; You're right, we of course need the negatives as well. Assuming the 2 coordinates are obtained through a rescaling x->2*x-1, if __sq is not exactly 0, it must be between 2^-103 and 1 (for ieee double), so I am not sure hypot gains that much (at least in my mind hypot was mostly a gain close to 0 or infinity, but maybe it has more advantages). It can only hurt speed though, so not a big issue. Depending on how similar in size the two values are, not using hypot() can drop quite a few bits. Especially with the scaling through division this error can be noticeable. Better be sure. Maybe at some point I have time to investigate the worst case scenario for the numbers in question but until this shows hypot isn't needed it's best to leave it in. I've committed the patch. Hi, This causes lots of tests under libstdc++-v3/testsuite/ext/ failed on aarch64-none-elf and arm-none-eabi with below log message. Please follow our usual rules, don't unconditionally use C99 functions like std::hypot, use _GLIBCXX_USE_C99_MATH_TR1 (many examples, eg in random.tcc). Thanks! Paolo.
Re: [PATCH] libstdc++: add uniform on sphere distribution
On Wed, Aug 13, 2014 at 3:36 PM, Bin.Cheng wrote:
> On Sun, Aug 10, 2014 at 2:00 AM, Ulrich Drepper wrote:
>> On Sat, Aug 9, 2014 at 1:40 PM, Marc Glisse wrote:
>>> __x = result_type(2.0) * __aurng() - 1.0;
>>
>> You're right, we of course need the negatives as well.
>>
>>> Assuming the 2 coordinates are obtained through a rescaling x->2*x-1, if
>>> __sq is not exactly 0, it must be between 2^-103 and 1 (for ieee
>>> double), so I am not sure hypot gains that much (at least in my mind
>>> hypot was mostly a gain close to 0 or infinity, but maybe it has more
>>> advantages). It can only hurt speed though, so not a big issue.
>>
>> Depending on how similar in size the two values are, not using hypot()
>> can drop quite a few bits. Especially with the scaling through
>> division this error can be noticeable. Better be sure. Maybe at some
>> point I have time to investigate the worst case scenario for the
>> numbers in question but until this shows hypot isn't needed it's best
>> to leave it in.
>>
>> I've committed the patch.
>
> Hi,
>
> This causes lots of tests under libstdc++-v3/testsuite/ext/ failed on
> aarch64-none-elf and arm-none-eabi with below log message.
>
> In file included from
> /home/binche01/work/oban-work/build-aarch64-none-elf/obj/gcc2/aarch64-none-elf/libstdc++-v3/include/ext/random:3494:0,
> from
> /home/binche01/work/oban-work/src/gcc/libstdc++-v3/testsuite/ext/random/arcsine_distribution/cons/default.cc:23:
> /home/binche01/work/oban-work/src/gcc/libstdc++-v3/include/ext/random.tcc:
> In member function
> '__gnu_cxx::{anonymous}::uniform_on_sphere_helper<2ul,
> _RealType>::result_type
> __gnu_cxx::{anonymous}::uniform_on_sphere_helper<2ul,
> _RealType>::operator()(_NormalDistribution&,
> _UniformRandomNumberGenerator&)':
> /home/binche01/work/oban-work/src/gcc/libstdc++-v3/include/ext/random.tcc:1609:18:
> error: 'hypot' is not a member of 'std'
> auto __norm = std::hypot(__ret[0], __ret[1]);
> ^
> /home/binche01/work/oban-work/src/gcc/libstdc++-v3/include/ext/random.tcc:1609:18:
> note: suggested alternative:
> In file included from
> /home/binche01/work/oban-work/build-aarch64-none-elf/obj/gcc2/aarch64-none-elf/libstdc++-v3/include/cmath:44:0,
> from
> /home/binche01/work/oban-work/build-aarch64-none-elf/obj/gcc2/aarch64-none-elf/libstdc++-v3/include/random:38,
> from
> /home/binche01/work/oban-work/build-aarch64-none-elf/obj/gcc2/aarch64-none-elf/libstdc++-v3/include/ext/random:38,
> from
> /home/binche01/work/oban-work/src/gcc/libstdc++-v3/testsuite/ext/random/arcsine_distribution/cons/default.cc:23:
> /home/binche01/work/oban-work/target-aarch64-none-elf/aarch64-none-elf/include/math.h:306:15:
> note: 'hypot'
> extern double hypot _PARAMS((double, double));
>^
>
> FAIL: ext/random/arcsine_distribution/cons/default.cc (test for excess errors)
> Excess errors:
> /home/binche01/work/oban-work/src/gcc/libstdc++-v3/include/ext/random.tcc:1609:18:
> error: 'hypot' is not a member of 'std'
>
https://gcc.gnu.org/bugzilla/show_bug.cgi?id=62118 filed.
Thanks,
bin
Re: [PATCH] libstdc++: add uniform on sphere distribution
On Sun, Aug 10, 2014 at 2:00 AM, Ulrich Drepper wrote:
> On Sat, Aug 9, 2014 at 1:40 PM, Marc Glisse wrote:
>> __x = result_type(2.0) * __aurng() - 1.0;
>
> You're right, we of course need the negatives as well.
>
>> Assuming the 2 coordinates are obtained through a rescaling x->2*x-1, if
>> __sq is not exactly 0, it must be between 2^-103 and 1 (for ieee
>> double), so I am not sure hypot gains that much (at least in my mind
>> hypot was mostly a gain close to 0 or infinity, but maybe it has more
>> advantages). It can only hurt speed though, so not a big issue.
>
> Depending on how similar in size the two values are, not using hypot()
> can drop quite a few bits. Especially with the scaling through
> division this error can be noticeable. Better be sure. Maybe at some
> point I have time to investigate the worst case scenario for the
> numbers in question but until this shows hypot isn't needed it's best
> to leave it in.
>
> I've committed the patch.
Hi,
This causes lots of tests under libstdc++-v3/testsuite/ext/ failed on
aarch64-none-elf and arm-none-eabi with below log message.
In file included from
/home/binche01/work/oban-work/build-aarch64-none-elf/obj/gcc2/aarch64-none-elf/libstdc++-v3/include/ext/random:3494:0,
from
/home/binche01/work/oban-work/src/gcc/libstdc++-v3/testsuite/ext/random/arcsine_distribution/cons/default.cc:23:
/home/binche01/work/oban-work/src/gcc/libstdc++-v3/include/ext/random.tcc:
In member function
'__gnu_cxx::{anonymous}::uniform_on_sphere_helper<2ul,
_RealType>::result_type
__gnu_cxx::{anonymous}::uniform_on_sphere_helper<2ul,
_RealType>::operator()(_NormalDistribution&,
_UniformRandomNumberGenerator&)':
/home/binche01/work/oban-work/src/gcc/libstdc++-v3/include/ext/random.tcc:1609:18:
error: 'hypot' is not a member of 'std'
auto __norm = std::hypot(__ret[0], __ret[1]);
^
/home/binche01/work/oban-work/src/gcc/libstdc++-v3/include/ext/random.tcc:1609:18:
note: suggested alternative:
In file included from
/home/binche01/work/oban-work/build-aarch64-none-elf/obj/gcc2/aarch64-none-elf/libstdc++-v3/include/cmath:44:0,
from
/home/binche01/work/oban-work/build-aarch64-none-elf/obj/gcc2/aarch64-none-elf/libstdc++-v3/include/random:38,
from
/home/binche01/work/oban-work/build-aarch64-none-elf/obj/gcc2/aarch64-none-elf/libstdc++-v3/include/ext/random:38,
from
/home/binche01/work/oban-work/src/gcc/libstdc++-v3/testsuite/ext/random/arcsine_distribution/cons/default.cc:23:
/home/binche01/work/oban-work/target-aarch64-none-elf/aarch64-none-elf/include/math.h:306:15:
note: 'hypot'
extern double hypot _PARAMS((double, double));
^
FAIL: ext/random/arcsine_distribution/cons/default.cc (test for excess errors)
Excess errors:
/home/binche01/work/oban-work/src/gcc/libstdc++-v3/include/ext/random.tcc:1609:18:
error: 'hypot' is not a member of 'std'
Thanks,
bin
Re: [PATCH] libstdc++: add uniform on sphere distribution
On Sat, Aug 9, 2014 at 1:40 PM, Marc Glisse wrote: > __x = result_type(2.0) * __aurng() - 1.0; You're right, we of course need the negatives as well. > Assuming the 2 coordinates are obtained through a rescaling x->2*x-1, if > __sq is not exactly 0, it must be between 2^-103 and 1 (for ieee > double), so I am not sure hypot gains that much (at least in my mind > hypot was mostly a gain close to 0 or infinity, but maybe it has more > advantages). It can only hurt speed though, so not a big issue. Depending on how similar in size the two values are, not using hypot() can drop quite a few bits. Especially with the scaling through division this error can be noticeable. Better be sure. Maybe at some point I have time to investigate the worst case scenario for the numbers in question but until this shows hypot isn't needed it's best to leave it in. I've committed the patch.
Re: [PATCH] libstdc++: add uniform on sphere distribution
On Sat, 9 Aug 2014, Ulrich Drepper wrote:
How about the patch below?
Looks good, with two details:
+template
+ class uniform_on_sphere_helper<2, _RealType>
+ {
+ typedef typename uniform_on_sphere_distribution<2, _RealType>::
+ result_type result_type;
+
+ public:
+ template
+ result_type operator()(_NormalDistribution&,
+ _UniformRandomNumberGenerator& __urng)
+{
+ result_type __ret;
+ _RealType __sq;
+ std::__detail::_Adaptor<_UniformRandomNumberGenerator,
+ _RealType> __aurng(__urng);
+
+ do
+ {
+ __ret[0] = __aurng();
I must be missing something obvious, but normal_distribution uses:
__x = result_type(2.0) * __aurng() - 1.0;
to get a number between -1 and 1, and I don't see where you do this
rescaling. Does __aurng() already return a number in the right interval in
this context? If so we may want to update our naming of variables to make
that clearer.
+ __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]);
Assuming the 2 coordinates are obtained through a rescaling x->2*x-1, if
__sq is not exactly 0, it must be between 2^-103 and 1 (for ieee
double), so I am not sure hypot gains that much (at least in my mind
hypot was mostly a gain close to 0 or infinity, but maybe it has more
advantages). It can only hurt speed though, so not a big issue.
--
Marc Glisse
Re: [PATCH] libstdc++: add uniform on sphere distribution
On 08/09/2014 11:33 AM, Marc Glisse wrote: On Sat, 9 Aug 2014, Ulrich Drepper wrote: If you are going to specialize for dim 2, I imagine you won't be computing normal distributions, you will only generate a point uniformy in a square and reject it if it is not in the ball? (interestingly enough this is used as a subroutine by the implementation of normal_distribution) We need to be *on* the circle, not inside. Yes, you still need the normalization step (divide by the norm). It works whether you start from a uniform distribution in the disk or from a Gaussian in the plane, but the first one is easier to generate (generate points in a square until you get one in the disk). When the dimension becomes higher, the probability that a point in the cube actually belongs to the ball decreases, and it becomes cheaper to use a Gaussian. If we pull in the implementation of normal you will just be able to use the two values that normal computes on the first, (and third, ...) calls without doing two calls. That and hypot would be a real win.
Re: [PATCH] libstdc++: add uniform on sphere distribution
Marc Glisse 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
* 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
+ class uniform_on_sphere_helper
+ {
+ typedef typename uniform_on_sphere_distribution<_Dimen,
_RealType>::result_type result_type;
+
+ public:
+ template
+ 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
+ class uniform_on_sphere_helper<2, _RealType>
+ {
+ typedef typename uniform_on_sphere_distribution<2, _RealType>::
+ result_type result_type;
+
+ public:
+ template
+ 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
template
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
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
//
Re: [PATCH] libstdc++: add uniform on sphere distribution
On Sat, 9 Aug 2014, Ulrich Drepper wrote: If you are going to specialize for dim 2, I imagine you won't be computing normal distributions, you will only generate a point uniformy in a square and reject it if it is not in the ball? (interestingly enough this is used as a subroutine by the implementation of normal_distribution) We need to be *on* the circle, not inside. Yes, you still need the normalization step (divide by the norm). It works whether you start from a uniform distribution in the disk or from a Gaussian in the plane, but the first one is easier to generate (generate points in a square until you get one in the disk). When the dimension becomes higher, the probability that a point in the cube actually belongs to the ball decreases, and it becomes cheaper to use a Gaussian. -- Marc Glisse
Re: [PATCH] libstdc++: add uniform on sphere distribution
On Sat, Aug 9, 2014 at 8:34 AM, Marc Glisse wrote: > Oh, a comment saying exactly what you just said would be fine with me (or > even nothing). We might at some point use a different method than Box-Muller sampling so I'm OK with the test. > If you are going to specialize for dim 2, I imagine you won't be computing > normal distributions, you will only generate a point uniformy in a square > and reject it if it is not in the ball? (interestingly enough this is used > as a subroutine by the implementation of normal_distribution) We need to be *on* the circle, not inside. We'll still have to follow the algorithm unless I miss something. With reasonable probability we cannot generate those numbers directly from a uniform source. What is optimized is just the norm computation.
Re: [PATCH] libstdc++: add uniform on sphere distribution
On Sat, 9 Aug 2014, Ulrich Drepper wrote: On Sat, Aug 9, 2014 at 3:15 AM, Marc Glisse wrote: While there, do we want to also reject infinite norms? I would have done: while (__sum < small || __sum > large) but testing exactly for 0 and infinity seems good enough. I guess the squaring can theoretically overflow and produce infinity. It will never happen with the way we generate normally distributed numbers, though. These values are always so unlikely that it is OK that the algorithms cannot return them. If you insist I'll add a test for infinity. Oh, a comment saying exactly what you just said would be fine with me (or even nothing). The other change (which would eliminate the necessity for this test in a special case) is to use hypot for _Dimen==2. This might be a case common enough to warrant that little bit of extra text. I'll prepare a patch. If you are going to specialize for dim 2, I imagine you won't be computing normal distributions, you will only generate a point uniformy in a square and reject it if it is not in the ball? (interestingly enough this is used as a subroutine by the implementation of normal_distribution) -- Marc Glisse
Re: [PATCH] libstdc++: add uniform on sphere distribution
On Sat, Aug 9, 2014 at 3:15 AM, Marc Glisse wrote: > While there, do we want to also reject infinite norms? > I would have done: while (__sum < small || __sum > large) > but testing exactly for 0 and infinity seems good enough. I guess the squaring can theoretically overflow and produce infinity. It will never happen with the way we generate normally distributed numbers, though. These values are always so unlikely that it is OK that the algorithms cannot return them. If you insist I'll add a test for infinity. The other change (which would eliminate the necessity for this test in a special case) is to use hypot for _Dimen==2. This might be a case common enough to warrant that little bit of extra text. I'll prepare a patch.
Re: [PATCH] libstdc++: add uniform on sphere distribution
On Fri, 8 Aug 2014, Ulrich Drepper wrote: Now also for a fix of the sphere distribution. Unless someone objects I'll check in the patch below. 2014-08-08 Ulrich Drepper * include/ext/random.tcc (uniform_on_sphere_distribution::__generate_impl): Reject vectors with norm zero. While there, do we want to also reject infinite norms? I would have done: while (__sum < small || __sum > large) but testing exactly for 0 and infinity seems good enough. -- Marc Glisse
Re: [PATCH] libstdc++: add uniform on sphere distribution
Jonathan Wakely writes:
> On 23/07/14 11:58 +0200, Marc Glisse wrote:
> As an aside, we already have divide-by-zero bugs in , it
> would be nice if someone could look at that.
>
> https://gcc.gnu.org/bugzilla/show_bug.cgi?id=60037
Sorry, it took a while to get back to tihs and now the referenced bug is
already fixed. Good. Now also for a fix of the sphere distribution.
Unless someone objects I'll check in the patch below.
2014-08-08 Ulrich Drepper
* include/ext/random.tcc
(uniform_on_sphere_distribution::__generate_impl): Reject
vectors with norm zero.
diff --git a/libstdc++-v3/include/ext/random.tcc
b/libstdc++-v3/include/ext/random.tcc
index 05361d8..d1f0b9c 100644
--- a/libstdc++-v3/include/ext/random.tcc
+++ b/libstdc++-v3/include/ext/random.tcc
@@ -1548,13 +1548,21 @@ _GLIBCXX_BEGIN_NAMESPACE_VERSION
const param_type& __p)
{
result_type __ret;
- _RealType __sum = _RealType(0);
+ _RealType __norm;
+
+ do
+ {
+ _RealType __sum = _RealType(0);
+
+ std::generate(__ret.begin(), __ret.end(),
+ [&__urng, &__sum, this](){
+ _RealType __t = _M_nd(__urng);
+ __sum += __t * __t;
+ return __t; });
+ __norm = std::sqrt(__sum);
+ }
+ while (__norm == _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; });
Re: [PATCH] libstdc++: add uniform on sphere distribution
On 24 July 2014 22:15, Ulrich Drepper wrote: > On Wed, Jul 23, 2014 at 6:29 AM, Jonathan Wakely wrote: >> As an aside, we already have divide-by-zero bugs in , it >> would be nice if someone could look at that. > > I'll take a look at this soon. That would be great, thanks!
Re: [PATCH] libstdc++: add uniform on sphere distribution
On Wed, Jul 23, 2014 at 6:29 AM, Jonathan Wakely wrote: > As an aside, we already have divide-by-zero bugs in , it > would be nice if someone could look at that. I'll take a look at this soon.
Re: [PATCH] libstdc++: add uniform on sphere distribution
On 23/07/14 11:58 +0200, Marc Glisse wrote: * can't we end up dividing by 0 if all values of the normal distribution happen to be 0? As an aside, we already have divide-by-zero bugs in , it would be nice if someone could look at that. https://gcc.gnu.org/bugzilla/show_bug.cgi?id=60037
Re: [PATCH] libstdc++: add uniform on sphere distribution
On Sat, 12 Jul 2014, Ulrich Drepper wrote: Ed's submission of the logistic regression distribution caused problems for me because, like Ed, I have changes to the header in my tree for a long time. Time to submit them. This first one is a new distribution. It generates coordinates for random points on a unit sphere in arbitrarily many dimensions. This distribution by itself is useful but if I get some other code fully implemented it will also form the basis for yet another, more sophisticated distribution. Hello, I have 2 questions : * can't we end up dividing by 0 if all values of the normal distribution happen to be 0? * should the implementation be specialized for small dimensions to avoid the normal distributions and instead generate points in a square/cube until they fall in the disk/ball? -- Marc Glisse
Re: [PATCH] libstdc++: add uniform on sphere distribution
On Wed, Jul 16, 2014 at 11:01 AM, Paolo Carlini wrote: > Right. And reset too. I'm going to test and apply the below. Sorry for not reacting to this, I was away from the machines I could have done that work on.
Re: [PATCH] libstdc++: add uniform on sphere distribution
Hi,
On 07/16/2014 01:05 PM, Ed Smith-Rowland wrote:
One thing we all forgot: the operator== is also non-trivial because it
needs to compare _M_n.
Right. And reset too. I'm going to test and apply the below.
Thanks,
Paolo.
///
2014-07-16 Paolo Carlini
* include/ext/random: Minor formatting and cosmetic tweaks.
(uniform_on_sphere_distribution<>::operator==
(const uniform_on_sphere_distribution&,
const uniform_on_sphere_distribution&)): Compare the _M_nds.
(uniform_on_sphere_distribution<>::reset): Reset _M_nd.
(operator!=(const uniform_on_sphere_distribution&,
const uniform_on_sphere_distribution&)): Adjust.
* include/ext/random.tcc: Minor cosmetc tweaks.
Index: include/ext/random
===
--- include/ext/random (revision 212581)
+++ include/ext/random (working copy)
@@ -598,7 +598,7 @@
inline bool
operator!=(const __gnu_cxx::beta_distribution<_RealType>& __d1,
const __gnu_cxx::beta_distribution<_RealType>& __d2)
- { return !(__d1 == __d2); }
+{ return !(__d1 == __d2); }
/**
@@ -2575,7 +2575,7 @@
inline bool
operator!=(const __gnu_cxx::triangular_distribution<_RealType>& __d1,
const __gnu_cxx::triangular_distribution<_RealType>& __d2)
- { return !(__d1 == __d2); }
+{ return !(__d1 == __d2); }
/**
@@ -2810,7 +2810,7 @@
inline bool
operator!=(const __gnu_cxx::von_mises_distribution<_RealType>& __d1,
const __gnu_cxx::von_mises_distribution<_RealType>& __d2)
- { return !(__d1 == __d2); }
+{ return !(__d1 == __d2); }
/**
@@ -3328,12 +3328,12 @@
*/
explicit
uniform_on_sphere_distribution()
- : _M_param(), _M_n(_RealType(0), _RealType(1))
+ : _M_param(), _M_nd()
{ }
explicit
uniform_on_sphere_distribution(const param_type& __p)
- : _M_param(__p), _M_n(_RealType(0), _RealType(1))
+ : _M_param(__p), _M_nd()
{ }
/**
@@ -3341,7 +3341,7 @@
*/
void
reset()
- { }
+ { _M_nd.reset(); }
/**
* @brief Returns the parameter set of the distribution.
@@ -3425,14 +3425,15 @@
friend bool
operator==(const uniform_on_sphere_distribution& __d1,
const uniform_on_sphere_distribution& __d2)
- { return true; }
+ { return __d1._M_nd == __d2._M_nd; }
/**
- * @brief Inserts a %uniform_on_sphere_distribution random number
distribution
- * @p __x into the output stream @p __os.
+ * @brief Inserts a %uniform_on_sphere_distribution random number
+ *distribution @p __x into the output stream @p __os.
*
* @param __os An output stream.
- * @param __x A %uniform_on_sphere_distribution random number
distribution.
+ * @param __x A %uniform_on_sphere_distribution random number
+ * distribution.
*
* @returns The output stream with the state of @p __x inserted or in
* an error state.
@@ -3446,11 +3447,13 @@
__x);
/**
- * @brief Extracts a %uniform_on_sphere_distribution random number
distribution
+ * @brief Extracts a %uniform_on_sphere_distribution random number
+ *distribution
* @p __x from the input stream @p __is.
*
* @param __is An input stream.
- * @param __x A %uniform_on_sphere_distribution random number generator
engine.
+ * @param __x A %uniform_on_sphere_distribution random number
+ * generator engine.
*
* @returns The input stream with @p __x extracted or in an error state.
*/
@@ -3470,7 +3473,7 @@
const param_type& __p);
param_type _M_param;
- std::normal_distribution<_RealType> _M_n;
+ std::normal_distribution<_RealType> _M_nd;
};
/**
@@ -3482,7 +3485,7 @@
_RealType>& __d1,
const __gnu_cxx::uniform_on_sphere_distribution<_Dimen,
_RealType>& __d2)
- { return false; }
+{ return !(__d1 == __d2); }
_GLIBCXX_END_NAMESPACE_VERSION
} // namespace __gnu_cxx
Index: include/ext/random.tcc
===
--- include/ext/random.tcc (revision 212581)
+++ include/ext/random.tcc (working copy)
@@ -1551,7 +1551,7 @@
_RealType __sum = _RealType(0);
std::generate(__ret.begin(), __ret.end(),
- [&__urng, &__sum, this](){ _RealType __t = _M_n(__urng);
+ [&__urng, &__sum, this](){ _RealType __t = _M_nd(__urng);
__sum += __t * __t;
return __t; });
auto __norm = std::sqrt(__sum);
@@ -1583,8 +1583,7 @@
const __gnu_cxx::uniform_on_sphere_distrib
Re: [PATCH] libstdc++: add uniform on sphere distribution
On 07/14/2014 04:14 AM, Paolo Carlini wrote: Hi, On 07/14/2014 09:58 AM, Andreas Schwab wrote: FAIL: ext/random/arcsine_distribution/cons/default.cc (test for excess errors) Excess errors: /daten/aranym/gcc/gcc-20140714/libstdc++-v3/include/ext/random.tcc:1587:22: error: '_M_n' was not declared in this scope /daten/aranym/gcc/gcc-20140714/libstdc++-v3/include/ext/random.tcc:1598:22: error: '_M_n' was not declared in this scope Ulrich please fix those _M_n to __x._M_n and remove the obsolete comments. We could also add the usual minimal set of tests exercising serialization. Thanks. Paolo. All, One thing we all forgot: the operator== is also non-trivial because it needs to compare _M_n. Operators must be guaranteed to produce the same numbers in future when fed with identical urngs if they compare equal. Ulrich, do you want to do this? There are also old comments in the inserters and extractors. Ed
Re: [PATCH] libstdc++: add uniform on sphere distribution
> I've checked in a patch to save the _M_n state. If it is r212496, did you test it? It breaks many tests (see https://gcc.gnu.org/ml/gcc-regression/2014-07/msg00213.html) /opt/gcc/work/libstdc++-v3/include/ext/random.tcc: In function 'std::basic_ostream<_CharT, _Traits>& __gnu_cxx::operator<<(std::basic_ostream<_CharT, _Traits>&, const __gnu_cxx::uniform_on_sphere_distribution<_Dimen, _RealType>&)': /opt/gcc/work/libstdc++-v3/include/ext/random.tcc:1587:22: error: '_M_n' was not declared in this scope return __os << _M_n; ^ /opt/gcc/work/libstdc++-v3/include/ext/random.tcc: In function 'std::basic_istream<_CharT, _Traits>& __gnu_cxx::operator>>(std::basic_istream<_CharT, _Traits>&, __gnu_cxx::uniform_on_sphere_distribution<_Dimen, _RealType>&)': /opt/gcc/work/libstdc++-v3/include/ext/random.tcc:1598:22: error: '_M_n' was not declared in this scope return __is >> _M_n ^ TIA Dominique
Re: [PATCH] libstdc++: add uniform on sphere distribution
Hi, On 07/14/2014 09:58 AM, Andreas Schwab wrote: FAIL: ext/random/arcsine_distribution/cons/default.cc (test for excess errors) Excess errors: /daten/aranym/gcc/gcc-20140714/libstdc++-v3/include/ext/random.tcc:1587:22: error: '_M_n' was not declared in this scope /daten/aranym/gcc/gcc-20140714/libstdc++-v3/include/ext/random.tcc:1598:22: error: '_M_n' was not declared in this scope Ulrich please fix those _M_n to __x._M_n and remove the obsolete comments. We could also add the usual minimal set of tests exercising serialization. Thanks. Paolo.
Re: [PATCH] libstdc++: add uniform on sphere distribution
FAIL: ext/random/arcsine_distribution/cons/default.cc (test for excess errors) Excess errors: /daten/aranym/gcc/gcc-20140714/libstdc++-v3/include/ext/random.tcc:1587:22: error: '_M_n' was not declared in this scope /daten/aranym/gcc/gcc-20140714/libstdc++-v3/include/ext/random.tcc:1598:22: error: '_M_n' was not declared in this scope Andreas. -- Andreas Schwab, SUSE Labs, sch...@suse.de GPG Key fingerprint = 0196 BAD8 1CE9 1970 F4BE 1748 E4D4 88E3 0EEA B9D7 "And now for something completely different."
Re: [PATCH] libstdc++: add uniform on sphere distribution
Hi, On 07/13/2014 06:44 PM, Ulrich Drepper wrote: But your 4th and 7th call example by itself is not a reason. Again, the input exclusively determined by the random numbers. Here, of course, the 4th and 7th use will produce different results. But this is not what the state of the distribution is supposed to capture. For that you'll have to save the state of the RNG as well. Yes, you are right. Saving the full state of the distribution solves only half of my hypothetical problem, but certainly you have to save it if you want to, say, reset the RNGs to a common state and get the same sequences of numbers after the os << x and is >> y pair. Paolo.
Re: [PATCH] libstdc++: add uniform on sphere distribution
On Sun, Jul 13, 2014 at 12:25 PM, Paolo Carlini wrote: > I don't think so. It depends on the past of the two different > uniform_on_sphere: each time each uniform_on_sphere calls _M_n(__urng) the > state of *its own* _M_n changes, evolves from the initial state. I indeed should use the normal_distribution operator<< and operator>> but I think for a different reason than you think. The way the normal_distribution is implemented produces two values at a time and saves the second for a latter call. So, yes, that implicit state has to be preserved and I should have followed what Ed said. But your 4th and 7th call example by itself is not a reason. Again, the input exclusively determined by the random numbers. Here, of course, the 4th and 7th use will produce different results. But this is not what the state of the distribution is supposed to capture. For that you'll have to save the state of the RNG as well. I've checked in a patch to save the _M_n state. Thanks.
Re: [PATCH] libstdc++: add uniform on sphere distribution
Hi, On 07/13/2014 06:18 PM, Ulrich Drepper wrote: On Sun, Jul 13, 2014 at 12:07 PM, Paolo Carlini wrote: Sorry, I still don't get it. When operator() of x and y, two uniform_on_sphere_distribution, call _M_n(__urng) and those _M_n have a different state, the numbers produced are in general different. Correct. But in the case of this distribution once you have the random numbers the remainder of the work is done by a fixed formula: v = (N(0,1), ..., N(0,1)) result = v / ||v||_2 That's it, nothing else to be done. If you have two calls of operator() of two different uniform_on_sphere objects and you pass to each an RNG object in exactly the same state you will get the same result. I don't think so. It depends on the past of the two different uniform_on_sphere: each time each uniform_on_sphere calls _M_n(__urng) the state of *its own* _M_n changes, evolves from the initial state. For instance, you call it 4 times on one and 7 times on another. The states of the two _M_n are different, and from that point the next numbers will be different. This is exactly what happens when you call os << x and is >> y, no object is destroyed, no object is constructed, and x has called _M_n(__urng) 4 times, and y has called it 7 times. Your inserters and extractors are trivial. The next numbers will be different. Note that, if I understand the standard correctly, we are comparing the numbers *from that point on*, not the whole story, whatever that may mean. In other terms, I'm saying the sequence of numbers produced after the 4th number isn't in general the same sequence of numbers produced after the 7th number. Paolo.
Re: [PATCH] libstdc++: add uniform on sphere distribution
On Sun, Jul 13, 2014 at 12:07 PM, Paolo Carlini wrote: > Sorry, I still don't get it. When operator() of x and y, two > uniform_on_sphere_distribution, call _M_n(__urng) and those _M_n have a > different state, the numbers produced are in general different. Correct. But in the case of this distribution once you have the random numbers the remainder of the work is done by a fixed formula: v = (N(0,1), ..., N(0,1)) result = v / ||v||_2 That's it, nothing else to be done. If you have two calls of operator() of two different uniform_on_sphere objects and you pass to each an RNG object in exactly the same state you will get the same result.
Re: [PATCH] libstdc++: add uniform on sphere distribution
Hi, On 07/13/2014 06:03 PM, Ulrich Drepper wrote: On Sun, Jul 13, 2014 at 11:43 AM, Paolo Carlini wrote: and I think: the normal distributions in x and y do have a non-trivial state (_M_saved, _M_saved_available) which, at any given moment, is different in x and y. Then the trivial inserter of x is called and the trivial extractor of y is called, nothing changes in y. I don't see how the following invocations of y(g) can produce the same sequence of numbers that would be produced by invocations of x(g). Remember: we are talking about distributions, not RNGs. The distribution has no parameters so given the same input (i.e., random byte sequences) it will create the same output all the time. Sorry, I still don't get it. When operator() of x and y, two uniform_on_sphere_distribution, call _M_n(__urng) and those _M_n have a different state, the numbers produced are in general different. For example, suppose _M_saved_available is true in both _M_n, everything is already "decided", no RNGs are involved. Paolo.
Re: [PATCH] libstdc++: add uniform on sphere distribution
On Sun, Jul 13, 2014 at 11:43 AM, Paolo Carlini wrote: > and I think: the normal distributions in x and y do have a non-trivial state > (_M_saved, _M_saved_available) which, at any given moment, is different in x > and y. Then the trivial inserter of x is called and the trivial extractor of > y is called, nothing changes in y. I don't see how the following invocations > of y(g) can produce the same sequence of numbers that would be produced by > invocations of x(g). Remember: we are talking about distributions, not RNGs. The distribution has no parameters so given the same input (i.e., random byte sequences) it will create the same output all the time.
Re: [PATCH] libstdc++: add uniform on sphere distribution
Hi, On 07/13/2014 04:11 PM, Ulrich Drepper wrote: On Sun, Jul 13, 2014 at 9:55 AM, Ed Smith-Rowland <3dw...@verizon.net> wrote: So I would just serialize _M_n here. It has fixed parameters. This would mean unnecessary work. When you try to use the parameter of the sphere distribution the normal distribution will be reset. So there really is no need here. At the cost of appearing a little dumb (and admittedly lately I'm not spending much time on , anyway), I'd like to make sure I *really* understand your reasoning... Thus I read this (26.5.1.6/6): "If a textual representation is written using os << x and that representation is restored into the same or a different object y of the same type using is >> y, repeated invocations of y(g) shall produce the same sequence of numbers as would repeated invocations of x(g)." and I think: the normal distributions in x and y do have a non-trivial state (_M_saved, _M_saved_available) which, at any given moment, is different in x and y. Then the trivial inserter of x is called and the trivial extractor of y is called, nothing changes in y. I don't see how the following invocations of y(g) can produce the same sequence of numbers that would be produced by invocations of x(g). Paolo.
Re: [PATCH] libstdc++: add uniform on sphere distribution
On 07/13/2014 10:11 AM, Ulrich Drepper wrote: On Sun, Jul 13, 2014 at 9:55 AM, Ed Smith-Rowland <3dw...@verizon.net> wrote: So I would just serialize _M_n here. It has fixed parameters. This would mean unnecessary work. When you try to use the parameter of the sphere distribution the normal distribution will be reset. So there really is no need here. The only problem would be if code couldn't handle the operators not writing/reading anything But I haven't seen anything like that. OK. I see it. Thanks.
Re: [PATCH] libstdc++: add uniform on sphere distribution
On Sun, Jul 13, 2014 at 9:55 AM, Ed Smith-Rowland <3dw...@verizon.net> wrote: > So I would just serialize _M_n here. It has fixed parameters. This would mean unnecessary work. When you try to use the parameter of the sphere distribution the normal distribution will be reset. So there really is no need here. The only problem would be if code couldn't handle the operators not writing/reading anything But I haven't seen anything like that.
Re: [PATCH] libstdc++: add uniform on sphere distribution
> are these dummy implementations intended? Yes. There is no state. The only parameter is the dimensionality which is a template parameter. We do often serialize underlying helper distributions, in your case the normal distribution _M_n. While the normal distribution mean and stddev are trivial for your case (not actually needing serialization) the normal distribution has these _M_saved, etc. shat we should store. So I would just serialize _M_n here. This is a great distribution. Thanks! Are you looking at the normal distribution equivalent of these?
Re: [PATCH] libstdc++: add uniform on sphere distribution
Hi, On 07/13/2014 12:04 PM, Ulrich Drepper wrote: On Sun, Jul 13, 2014 at 5:24 AM, Paolo Carlini wrote: are these dummy implementations intended? Yes. There is no state. The only parameter is the dimensionality which is a template parameter. Ah, interesting, didn't look close enough, sorry. Maybe add a one line comment in the bodies? Patch is Ok. Thanks! Paolo.
Re: [PATCH] libstdc++: add uniform on sphere distribution
On Sun, Jul 13, 2014 at 5:24 AM, Paolo Carlini wrote: > are these dummy implementations intended? Yes. There is no state. The only parameter is the dimensionality which is a template parameter.
Re: [PATCH] libstdc++: add uniform on sphere distribution
Hi,
On 07/13/2014 04:00 AM, Ulrich Drepper wrote:
+ template
+std::basic_ostream<_CharT, _Traits>&
+operator<<(std::basic_ostream<_CharT, _Traits>& __os,
+ const __gnu_cxx::uniform_on_sphere_distribution<_Dimen,
+ _RealType>& __x)
+{
+ return __os;
+}
+
+ template
+std::basic_istream<_CharT, _Traits>&
+operator>>(std::basic_istream<_CharT, _Traits>& __is,
+ __gnu_cxx::uniform_on_sphere_distribution<_Dimen,
+_RealType>& __x)
+{
+ return __is;
+}
are these dummy implementations intended?
Thanks,
Paolo.
Re: [PATCH] libstdc++: add uniform on sphere distribution
On Sat, Jul 12, 2014 at 10:00 PM, Ulrich Drepper wrote: > The patch is tested against the current tree without causing additional > problems. > > OK? Ignore the values.cc test case, it's a left-over from copying existing tests and modifying them for a new distribution. This test does not apply to the new distribution and therefore the entire file is gone.
