Mercurial > hg > octave-lyh
changeset 17266:213ee68b59da
Handle out-of-range values consistently when initializing random number generator.
* oct-rand.cc: New function double2uint32(). It normalizes a finite double by
taking it modulo 2^32, such that the result is in [0, 2^32). Nonfinites give
zero as result. This prevents compiler-dependent casting of a double whose
truncated value cannot be represented in a uint32_t.
* rand.cc: Add %!test block to verify behavior.
author | Philipp Kutin <philipp.kutin@gmail.com> |
---|---|
date | Thu, 15 Aug 2013 20:14:30 +0200 |
parents | ce55a8575874 |
children | a654c1eccd29 |
files | libinterp/corefcn/rand.cc liboctave/numeric/oct-rand.cc |
diffstat | 2 files changed, 38 insertions(+), 1 deletions(-) [+] |
line wrap: on
line diff
--- a/libinterp/corefcn/rand.cc +++ b/libinterp/corefcn/rand.cc @@ -533,6 +533,19 @@ %! endif */ +/* +%!# Test out-of-range values as rand() seeds. See oct-rand.cc: double2uint32(). +%!function v = __rand_sample__ (initval) +%! rand ("state", initval); +%! v = rand (1, 6); +%!endfunction +%! +%!assert (__rand_sample__ (0), __rand_sample__ (2^32)) +%!assert (__rand_sample__ (-2), __rand_sample__ (2^32-2)) +%!assert (__rand_sample__ (Inf), __rand_sample__ (NaN)) +%!assert (! isequal (__rand_sample__ (-1), __rand_sample__ (-2))) +*/ + static std::string current_distribution = octave_rand::distribution (); DEFUN (randn, args, ,
--- a/liboctave/numeric/oct-rand.cc +++ b/liboctave/numeric/oct-rand.cc @@ -663,6 +663,30 @@ return retval; } +// Guarantee reproducible conversion of negative initialization values to +// random number algorithm. Note that Matlab employs slightly different rules. +// 1) Seed saturates at 2^32-1 for any value larger than that. +// 2) NaN, Inf are translated to 2^32-1. +// 3) -Inf is translated to 0. +static uint32_t +double2uint32 (double d) +{ + uint32_t u; + static const double TWOUP32 = std::numeric_limits<uint32_t>::max() + 1.0; + + if (! xfinite (d)) + u = 0; + else + { + d = fmod (d, TWOUP32); + if (d < 0) + d += TWOUP32; + u = static_cast<uint32_t> (d); + } + + return u; +} + void octave_rand::set_internal_state (const ColumnVector& s) { @@ -672,7 +696,7 @@ OCTAVE_LOCAL_BUFFER (uint32_t, tmp, MT_N + 1); for (octave_idx_type i = 0; i < n; i++) - tmp[i] = static_cast<uint32_t> (s.elem (i)); + tmp[i] = double2uint32 (s.elem (i)); if (len == MT_N + 1 && tmp[MT_N] <= MT_N && tmp[MT_N] > 0) oct_set_state (tmp);