Back to home page

EIC code displayed by LXR

 
 

    


File indexing completed on 2025-01-18 09:51:13

0001 /* boost random/uniform_int_distribution.hpp header file
0002  *
0003  * Copyright Jens Maurer 2000-2001
0004  * Copyright Steven Watanabe 2011
0005  * Distributed under the Boost Software License, Version 1.0. (See
0006  * accompanying file LICENSE_1_0.txt or copy at
0007  * http://www.boost.org/LICENSE_1_0.txt)
0008  *
0009  * See http://www.boost.org for most recent version including documentation.
0010  *
0011  * $Id$
0012  *
0013  * Revision history
0014  *  2001-04-08  added min<max assertion (N. Becker)
0015  *  2001-02-18  moved to individual header files
0016  */
0017 
0018 #ifndef BOOST_RANDOM_UNIFORM_INT_DISTRIBUTION_HPP
0019 #define BOOST_RANDOM_UNIFORM_INT_DISTRIBUTION_HPP
0020 
0021 #include <iosfwd>
0022 #include <ios>
0023 #include <istream>
0024 #include <boost/config.hpp>
0025 #include <boost/limits.hpp>
0026 #include <boost/assert.hpp>
0027 #include <boost/random/detail/config.hpp>
0028 #include <boost/random/detail/operators.hpp>
0029 #include <boost/random/detail/uniform_int_float.hpp>
0030 #include <boost/random/detail/signed_unsigned_tools.hpp>
0031 #include <boost/random/traits.hpp>
0032 #include <boost/type_traits/integral_constant.hpp>
0033 #ifdef BOOST_NO_CXX11_EXPLICIT_CONVERSION_OPERATORS
0034 #include <boost/type_traits/conditional.hpp>
0035 #endif
0036 
0037 namespace boost {
0038 namespace random {
0039 namespace detail {
0040     
0041 
0042 #ifdef BOOST_MSVC
0043 #pragma warning(push)
0044 // disable division by zero warning, since we can't
0045 // actually divide by zero.
0046 #pragma warning(disable:4723)
0047 #endif
0048 
0049 template<class Engine, class T>
0050 T generate_uniform_int(
0051     Engine& eng, T min_value, T max_value,
0052     boost::true_type /** is_integral<Engine::result_type> */)
0053 {
0054     typedef T result_type;
0055     typedef typename boost::random::traits::make_unsigned_or_unbounded<T>::type range_type;
0056     typedef typename Engine::result_type base_result;
0057     // ranges are always unsigned or unbounded
0058     typedef typename boost::random::traits::make_unsigned_or_unbounded<base_result>::type base_unsigned;
0059     const range_type range = random::detail::subtract<result_type>()(max_value, min_value);
0060     const base_result bmin = (eng.min)();
0061     const base_unsigned brange =
0062       random::detail::subtract<base_result>()((eng.max)(), (eng.min)());
0063 
0064     if(range == 0) {
0065       return min_value;    
0066     } else if(brange == range) {
0067       // this will probably never happen in real life
0068       // basically nothing to do; just take care we don't overflow / underflow
0069       base_unsigned v = random::detail::subtract<base_result>()(eng(), bmin);
0070       return random::detail::add<base_unsigned, result_type>()(v, min_value);
0071     } else if(brange < range) {
0072       // use rejection method to handle things like 0..3 --> 0..4
0073       for(;;) {
0074         // concatenate several invocations of the base RNG
0075         // take extra care to avoid overflows
0076 
0077         //  limit == floor((range+1)/(brange+1))
0078         //  Therefore limit*(brange+1) <= range+1
0079         range_type limit;
0080         if(range == (std::numeric_limits<range_type>::max)()) {
0081           limit = range/(range_type(brange)+1);
0082           if(range % (range_type(brange)+1) == range_type(brange))
0083             ++limit;
0084         } else {
0085           limit = (range+1)/(range_type(brange)+1);
0086         }
0087 
0088         // We consider "result" as expressed to base (brange+1):
0089         // For every power of (brange+1), we determine a random factor
0090         range_type result = range_type(0);
0091         range_type mult = range_type(1);
0092 
0093         // loop invariants:
0094         //  result < mult
0095         //  mult <= range
0096         while(mult <= limit) {
0097           // Postcondition: result <= range, thus no overflow
0098           //
0099           // limit*(brange+1)<=range+1                   def. of limit       (1)
0100           // eng()-bmin<=brange                          eng() post.         (2)
0101           // and mult<=limit.                            loop condition      (3)
0102           // Therefore mult*(eng()-bmin+1)<=range+1      by (1),(2),(3)      (4)
0103           // Therefore mult*(eng()-bmin)+mult<=range+1   rearranging (4)     (5)
0104           // result<mult                                 loop invariant      (6)
0105           // Therefore result+mult*(eng()-bmin)<range+1  by (5), (6)         (7)
0106           //
0107           // Postcondition: result < mult*(brange+1)
0108           //
0109           // result<mult                                 loop invariant      (1)
0110           // eng()-bmin<=brange                          eng() post.         (2)
0111           // Therefore result+mult*(eng()-bmin) <
0112           //           mult+mult*(eng()-bmin)            by (1)              (3)
0113           // Therefore result+(eng()-bmin)*mult <
0114           //           mult+mult*brange                  by (2), (3)         (4)
0115           // Therefore result+(eng()-bmin)*mult <
0116           //           mult*(brange+1)                   by (4)
0117           result += static_cast<range_type>(static_cast<range_type>(random::detail::subtract<base_result>()(eng(), bmin)) * mult);
0118 
0119           // equivalent to (mult * (brange+1)) == range+1, but avoids overflow.
0120           if(mult * range_type(brange) == range - mult + 1) {
0121               // The destination range is an integer power of
0122               // the generator's range.
0123               return(result);
0124           }
0125 
0126           // Postcondition: mult <= range
0127           // 
0128           // limit*(brange+1)<=range+1                   def. of limit       (1)
0129           // mult<=limit                                 loop condition      (2)
0130           // Therefore mult*(brange+1)<=range+1          by (1), (2)         (3)
0131           // mult*(brange+1)!=range+1                    preceding if        (4)
0132           // Therefore mult*(brange+1)<range+1           by (3), (4)         (5)
0133           // 
0134           // Postcondition: result < mult
0135           //
0136           // See the second postcondition on the change to result. 
0137           mult *= range_type(brange)+range_type(1);
0138         }
0139         // loop postcondition: range/mult < brange+1
0140         //
0141         // mult > limit                                  loop condition      (1)
0142         // Suppose range/mult >= brange+1                Assumption          (2)
0143         // range >= mult*(brange+1)                      by (2)              (3)
0144         // range+1 > mult*(brange+1)                     by (3)              (4)
0145         // range+1 > (limit+1)*(brange+1)                by (1), (4)         (5)
0146         // (range+1)/(brange+1) > limit+1                by (5)              (6)
0147         // limit < floor((range+1)/(brange+1))           by (6)              (7)
0148         // limit==floor((range+1)/(brange+1))            def. of limit       (8)
0149         // not (2)                                       reductio            (9)
0150         //
0151         // loop postcondition: (range/mult)*mult+(mult-1) >= range
0152         //
0153         // (range/mult)*mult + range%mult == range       identity            (1)
0154         // range%mult < mult                             def. of %           (2)
0155         // (range/mult)*mult+mult > range                by (1), (2)         (3)
0156         // (range/mult)*mult+(mult-1) >= range           by (3)              (4)
0157         //
0158         // Note that the maximum value of result at this point is (mult-1),
0159         // so after this final step, we generate numbers that can be
0160         // at least as large as range.  We have to really careful to avoid
0161         // overflow in this final addition and in the rejection.  Anything
0162         // that overflows is larger than range and can thus be rejected.
0163 
0164         // range/mult < brange+1  -> no endless loop
0165         range_type result_increment =
0166             generate_uniform_int(
0167                 eng,
0168                 static_cast<range_type>(0),
0169                 static_cast<range_type>(range/mult),
0170                 boost::true_type());
0171         if(std::numeric_limits<range_type>::is_bounded && ((std::numeric_limits<range_type>::max)() / mult < result_increment)) {
0172           // The multiplcation would overflow.  Reject immediately.
0173           continue;
0174         }
0175         result_increment *= mult;
0176         // unsigned integers are guaranteed to wrap on overflow.
0177         result += result_increment;
0178         if(result < result_increment) {
0179           // The addition overflowed.  Reject.
0180           continue;
0181         }
0182         if(result > range) {
0183           // Too big.  Reject.
0184           continue;
0185         }
0186         return random::detail::add<range_type, result_type>()(result, min_value);
0187       }
0188     } else {                   // brange > range
0189 #ifdef BOOST_NO_CXX11_EXPLICIT_CONVERSION_OPERATORS
0190       typedef typename conditional<
0191          std::numeric_limits<range_type>::is_specialized && std::numeric_limits<base_unsigned>::is_specialized
0192          && (std::numeric_limits<range_type>::digits >= std::numeric_limits<base_unsigned>::digits),
0193          range_type, base_unsigned>::type mixed_range_type;
0194 #else
0195       typedef base_unsigned mixed_range_type;
0196 #endif
0197 
0198       mixed_range_type bucket_size;
0199       // it's safe to add 1 to range, as long as we cast it first,
0200       // because we know that it is less than brange.  However,
0201       // we do need to be careful not to cause overflow by adding 1
0202       // to brange.  We use mixed_range_type throughout for mixed
0203       // arithmetic between base_unsigned and range_type - in the case
0204       // that range_type has more bits than base_unsigned it is always
0205       // safe to use range_type for this albeit it may be more effient
0206       // to use base_unsigned.  The latter is a narrowing conversion though
0207       // which may be disallowed if range_type is a multiprecision type
0208       // and there are no explicit converison operators.
0209 
0210       if(brange == (std::numeric_limits<base_unsigned>::max)()) {
0211         bucket_size = static_cast<mixed_range_type>(brange) / (static_cast<mixed_range_type>(range)+1);
0212         if(static_cast<mixed_range_type>(brange) % (static_cast<mixed_range_type>(range)+1) == static_cast<mixed_range_type>(range)) {
0213           ++bucket_size;
0214         }
0215       } else {
0216         bucket_size = static_cast<mixed_range_type>(brange + 1) / (static_cast<mixed_range_type>(range)+1);
0217       }
0218       for(;;) {
0219         mixed_range_type result =
0220           random::detail::subtract<base_result>()(eng(), bmin);
0221         result /= bucket_size;
0222         // result and range are non-negative, and result is possibly larger
0223         // than range, so the cast is safe
0224         if(result <= static_cast<mixed_range_type>(range))
0225           return random::detail::add<mixed_range_type, result_type>()(result, min_value);
0226       }
0227     }
0228 }
0229 
0230 #ifdef BOOST_MSVC
0231 #pragma warning(pop)
0232 #endif
0233 
0234 template<class Engine, class T>
0235 inline T generate_uniform_int(
0236     Engine& eng, T min_value, T max_value,
0237     boost::false_type /** is_integral<Engine::result_type> */)
0238 {
0239     uniform_int_float<Engine> wrapper(eng);
0240     return generate_uniform_int(wrapper, min_value, max_value, boost::true_type());
0241 }
0242 
0243 template<class Engine, class T>
0244 inline T generate_uniform_int(Engine& eng, T min_value, T max_value)
0245 {
0246     typedef typename Engine::result_type base_result;
0247     return generate_uniform_int(eng, min_value, max_value,
0248         boost::random::traits::is_integral<base_result>());
0249 }
0250 
0251 }
0252 
0253 /**
0254  * The class template uniform_int_distribution models a \random_distribution.
0255  * On each invocation, it returns a random integer value uniformly
0256  * distributed in the set of integers {min, min+1, min+2, ..., max}.
0257  *
0258  * The template parameter IntType shall denote an integer-like value type.
0259  */
0260 template<class IntType = int>
0261 class uniform_int_distribution
0262 {
0263 public:
0264     typedef IntType input_type;
0265     typedef IntType result_type;
0266 
0267     class param_type
0268     {
0269     public:
0270 
0271         typedef uniform_int_distribution distribution_type;
0272 
0273         /**
0274          * Constructs the parameters of a uniform_int_distribution.
0275          *
0276          * Requires min <= max
0277          */
0278         explicit param_type(
0279             IntType min_arg = 0,
0280             IntType max_arg = (std::numeric_limits<IntType>::max)())
0281           : _min(min_arg), _max(max_arg)
0282         {
0283             BOOST_ASSERT(_min <= _max);
0284         }
0285 
0286         /** Returns the minimum value of the distribution. */
0287         IntType a() const { return _min; }
0288         /** Returns the maximum value of the distribution. */
0289         IntType b() const { return _max; }
0290 
0291         /** Writes the parameters to a @c std::ostream. */
0292         BOOST_RANDOM_DETAIL_OSTREAM_OPERATOR(os, param_type, parm)
0293         {
0294             os << parm._min << " " << parm._max;
0295             return os;
0296         }
0297 
0298         /** Reads the parameters from a @c std::istream. */
0299         BOOST_RANDOM_DETAIL_ISTREAM_OPERATOR(is, param_type, parm)
0300         {
0301             IntType min_in, max_in;
0302             if(is >> min_in >> std::ws >> max_in) {
0303                 if(min_in <= max_in) {
0304                     parm._min = min_in;
0305                     parm._max = max_in;
0306                 } else {
0307                     is.setstate(std::ios_base::failbit);
0308                 }
0309             }
0310             return is;
0311         }
0312 
0313         /** Returns true if the two sets of parameters are equal. */
0314         BOOST_RANDOM_DETAIL_EQUALITY_OPERATOR(param_type, lhs, rhs)
0315         { return lhs._min == rhs._min && lhs._max == rhs._max; }
0316 
0317         /** Returns true if the two sets of parameters are different. */
0318         BOOST_RANDOM_DETAIL_INEQUALITY_OPERATOR(param_type)
0319 
0320     private:
0321 
0322         IntType _min;
0323         IntType _max;
0324     };
0325 
0326     /**
0327      * Constructs a uniform_int_distribution. @c min and @c max are
0328      * the parameters of the distribution.
0329      *
0330      * Requires: min <= max
0331      */
0332     explicit uniform_int_distribution(
0333         IntType min_arg = 0,
0334         IntType max_arg = (std::numeric_limits<IntType>::max)())
0335       : _min(min_arg), _max(max_arg)
0336     {
0337         BOOST_ASSERT(min_arg <= max_arg);
0338     }
0339     /** Constructs a uniform_int_distribution from its parameters. */
0340     explicit uniform_int_distribution(const param_type& parm)
0341       : _min(parm.a()), _max(parm.b()) {}
0342 
0343     /**  Returns the minimum value of the distribution */
0344     IntType min BOOST_PREVENT_MACRO_SUBSTITUTION () const { return _min; }
0345     /**  Returns the maximum value of the distribution */
0346     IntType max BOOST_PREVENT_MACRO_SUBSTITUTION () const { return _max; }
0347 
0348     /**  Returns the minimum value of the distribution */
0349     IntType a() const { return _min; }
0350     /**  Returns the maximum value of the distribution */
0351     IntType b() const { return _max; }
0352 
0353     /** Returns the parameters of the distribution. */
0354     param_type param() const { return param_type(_min, _max); }
0355     /** Sets the parameters of the distribution. */
0356     void param(const param_type& parm)
0357     {
0358         _min = parm.a();
0359         _max = parm.b();
0360     }
0361 
0362     /**
0363      * Effects: Subsequent uses of the distribution do not depend
0364      * on values produced by any engine prior to invoking reset.
0365      */
0366     void reset() { }
0367 
0368     /** Returns an integer uniformly distributed in the range [min, max]. */
0369     template<class Engine>
0370     result_type operator()(Engine& eng) const
0371     { return detail::generate_uniform_int(eng, _min, _max); }
0372 
0373     /**
0374      * Returns an integer uniformly distributed in the range
0375      * [param.a(), param.b()].
0376      */
0377     template<class Engine>
0378     result_type operator()(Engine& eng, const param_type& parm) const
0379     { return detail::generate_uniform_int(eng, parm.a(), parm.b()); }
0380 
0381     /** Writes the distribution to a @c std::ostream. */
0382     BOOST_RANDOM_DETAIL_OSTREAM_OPERATOR(os, uniform_int_distribution, ud)
0383     {
0384         os << ud.param();
0385         return os;
0386     }
0387 
0388     /** Reads the distribution from a @c std::istream. */
0389     BOOST_RANDOM_DETAIL_ISTREAM_OPERATOR(is, uniform_int_distribution, ud)
0390     {
0391         param_type parm;
0392         if(is >> parm) {
0393             ud.param(parm);
0394         }
0395         return is;
0396     }
0397 
0398     /**
0399      * Returns true if the two distributions will produce identical sequences
0400      * of values given equal generators.
0401      */
0402     BOOST_RANDOM_DETAIL_EQUALITY_OPERATOR(uniform_int_distribution, lhs, rhs)
0403     { return lhs._min == rhs._min && lhs._max == rhs._max; }
0404     
0405     /**
0406      * Returns true if the two distributions may produce different sequences
0407      * of values given equal generators.
0408      */
0409     BOOST_RANDOM_DETAIL_INEQUALITY_OPERATOR(uniform_int_distribution)
0410 
0411 private:
0412     IntType _min;
0413     IntType _max;
0414 };
0415 
0416 } // namespace random
0417 } // namespace boost
0418 
0419 #endif // BOOST_RANDOM_UNIFORM_INT_HPP