bug #1584: Improve random (avoid undefined behavior).

This commit is contained in:
Alexey Frunze 2018-08-08 20:19:32 -07:00
parent 5747288676
commit 050bcf6126
2 changed files with 40 additions and 12 deletions

View File

@ -684,20 +684,27 @@ struct random_default_impl<Scalar, false, true>
{
static inline Scalar run(const Scalar& x, const Scalar& y)
{
typedef typename conditional<NumTraits<Scalar>::IsSigned,std::ptrdiff_t,std::size_t>::type ScalarX;
if(y<x)
if (y <= x)
return x;
// the following difference might overflow on a 32 bits system,
// but since y>=x the result converted to an unsigned long is still correct.
std::size_t range = ScalarX(y)-ScalarX(x);
std::size_t offset = 0;
// rejection sampling
std::size_t divisor = 1;
std::size_t multiplier = 1;
if(range<RAND_MAX) divisor = (std::size_t(RAND_MAX)+1)/(range+1);
else multiplier = 1 + range/(std::size_t(RAND_MAX)+1);
// ScalarU is the unsigned counterpart of Scalar, possibly Scalar itself.
typedef typename make_unsigned<Scalar>::type ScalarU;
// ScalarX is the widest of ScalarU and unsigned int.
// We'll deal only with ScalarX and unsigned int below thus avoiding signed
// types and arithmetic and signed overflows (which are undefined behavior).
typedef typename conditional<(ScalarU(-1) > unsigned(-1)), ScalarU, unsigned>::type ScalarX;
// The following difference doesn't overflow, provided our integer types are two's
// complement and have the same number of padding bits in signed and unsigned variants.
// This is the case in most modern implementations of C++.
ScalarX range = ScalarX(y) - ScalarX(x);
ScalarX offset = 0;
ScalarX divisor = 1;
ScalarX multiplier = 1;
const unsigned rand_max = RAND_MAX;
if (range <= rand_max) divisor = (rand_max + 1) / (range + 1);
else multiplier = 1 + range / (rand_max + 1);
// Rejection sampling.
do {
offset = (std::size_t(std::rand()) * multiplier) / divisor;
offset = (unsigned(std::rand()) * multiplier) / divisor;
} while (offset > range);
return Scalar(ScalarX(x) + offset);
}

View File

@ -128,6 +128,27 @@ template<> struct is_integral<unsigned __int64> { enum { value = true }; }
#endif
#endif
#if EIGEN_HAS_CXX11
using std::make_unsigned;
#else
// TODO: Possibly improve this implementation of make_unsigned.
// It is currently used only by
// template<typename Scalar> struct random_default_impl<Scalar, false, true>.
template<typename> struct make_unsigned;
template<> struct make_unsigned<char> { typedef unsigned char type; };
template<> struct make_unsigned<signed char> { typedef unsigned char type; };
template<> struct make_unsigned<unsigned char> { typedef unsigned char type; };
template<> struct make_unsigned<signed short> { typedef unsigned short type; };
template<> struct make_unsigned<unsigned short> { typedef unsigned short type; };
template<> struct make_unsigned<signed int> { typedef unsigned int type; };
template<> struct make_unsigned<unsigned int> { typedef unsigned int type; };
template<> struct make_unsigned<signed long> { typedef unsigned long type; };
template<> struct make_unsigned<unsigned long> { typedef unsigned long type; };
#if EIGEN_COMP_MSVC
template<> struct make_unsigned<signed __int64> { typedef unsigned __int64 type; };
template<> struct make_unsigned<unsigned __int64> { typedef unsigned __int64 type; };
#endif
#endif
template <typename T> struct add_const { typedef const T type; };
template <typename T> struct add_const<T&> { typedef T& type; };