Chris@16: /* boost random/subtract_with_carry.hpp header file
Chris@16:  *
Chris@16:  * Copyright Jens Maurer 2002
Chris@16:  * Distributed under the Boost Software License, Version 1.0. (See
Chris@16:  * accompanying file LICENSE_1_0.txt or copy at
Chris@16:  * http://www.boost.org/LICENSE_1_0.txt)
Chris@16:  *
Chris@16:  * See http://www.boost.org for most recent version including documentation.
Chris@16:  *
Chris@101:  * $Id$
Chris@16:  *
Chris@16:  * Revision history
Chris@16:  *  2002-03-02  created
Chris@16:  */
Chris@16: 
Chris@16: #ifndef BOOST_RANDOM_SUBTRACT_WITH_CARRY_HPP
Chris@16: #define BOOST_RANDOM_SUBTRACT_WITH_CARRY_HPP
Chris@16: 
Chris@16: #include <boost/config/no_tr1/cmath.hpp>         // std::pow
Chris@16: #include <iostream>
Chris@16: #include <algorithm>     // std::equal
Chris@16: #include <stdexcept>
Chris@16: #include <boost/config.hpp>
Chris@16: #include <boost/limits.hpp>
Chris@16: #include <boost/cstdint.hpp>
Chris@16: #include <boost/static_assert.hpp>
Chris@16: #include <boost/integer/static_log2.hpp>
Chris@16: #include <boost/integer/integer_mask.hpp>
Chris@16: #include <boost/detail/workaround.hpp>
Chris@16: #include <boost/random/detail/config.hpp>
Chris@16: #include <boost/random/detail/seed.hpp>
Chris@16: #include <boost/random/detail/operators.hpp>
Chris@16: #include <boost/random/detail/seed_impl.hpp>
Chris@16: #include <boost/random/detail/generator_seed_seq.hpp>
Chris@16: #include <boost/random/linear_congruential.hpp>
Chris@16: 
Chris@16: 
Chris@16: namespace boost {
Chris@16: namespace random {
Chris@16: 
Chris@16: namespace detail {
Chris@16:    
Chris@16: struct subtract_with_carry_discard
Chris@16: {
Chris@16:     template<class Engine>
Chris@16:     static void apply(Engine& eng, boost::uintmax_t z)
Chris@16:     {
Chris@16:         typedef typename Engine::result_type IntType;
Chris@16:         const std::size_t short_lag = Engine::short_lag;
Chris@16:         const std::size_t long_lag = Engine::long_lag;
Chris@16:         std::size_t k = eng.k;
Chris@16:         IntType carry = eng.carry;
Chris@16:         if(k != 0) {
Chris@16:             // increment k until it becomes 0.
Chris@16:             if(k < short_lag) {
Chris@16:                 std::size_t limit = (short_lag - k) < z?
Chris@16:                     short_lag : (k + static_cast<std::size_t>(z));
Chris@16:                 for(std::size_t j = k; j < limit; ++j) {
Chris@16:                     carry = eng.do_update(j, j + long_lag - short_lag, carry);
Chris@16:                 }
Chris@16:             }
Chris@16:             std::size_t limit = (long_lag - k) < z?
Chris@16:                 long_lag : (k + static_cast<std::size_t>(z));
Chris@16:             std::size_t start = (k < short_lag ? short_lag : k);
Chris@16:             for(std::size_t j = start; j < limit; ++j) {
Chris@16:                 carry = eng.do_update(j, j - short_lag, carry);
Chris@16:             }
Chris@16:         }
Chris@16: 
Chris@16:         k = ((z % long_lag) + k) % long_lag;
Chris@16: 
Chris@16:         if(k < z) {
Chris@16:             // main loop: update full blocks from k = 0 to long_lag
Chris@16:             for(std::size_t i = 0; i < (z - k) / long_lag; ++i) {
Chris@16:                 for(std::size_t j = 0; j < short_lag; ++j) {
Chris@16:                     carry = eng.do_update(j, j + long_lag - short_lag, carry);
Chris@16:                 }
Chris@16:                 for(std::size_t j = short_lag; j < long_lag; ++j) {
Chris@16:                     carry = eng.do_update(j, j - short_lag, carry);
Chris@16:                 }
Chris@16:             }
Chris@16: 
Chris@16:             // Update the last partial block
Chris@16:             std::size_t limit = short_lag < k? short_lag : k; 
Chris@16:             for(std::size_t j = 0; j < limit; ++j) {
Chris@16:                 carry = eng.do_update(j, j + long_lag - short_lag, carry);
Chris@16:             }
Chris@16:             for(std::size_t j = short_lag; j < k; ++j) {
Chris@16:                 carry = eng.do_update(j, j - short_lag, carry);
Chris@16:             }
Chris@16:         }
Chris@16:         eng.carry = carry;
Chris@16:         eng.k = k;
Chris@16:     }
Chris@16: };
Chris@16: 
Chris@16: }
Chris@16: 
Chris@16: /**
Chris@16:  * Instantiations of @c subtract_with_carry_engine model a
Chris@16:  * \pseudo_random_number_generator.  The algorithm is
Chris@16:  * described in
Chris@16:  *
Chris@16:  *  @blockquote
Chris@16:  *  "A New Class of Random Number Generators", George
Chris@16:  *  Marsaglia and Arif Zaman, Annals of Applied Probability,
Chris@16:  *  Volume 1, Number 3 (1991), 462-480.
Chris@16:  *  @endblockquote
Chris@16:  */
Chris@16: template<class IntType, std::size_t w, std::size_t s, std::size_t r>
Chris@16: class subtract_with_carry_engine
Chris@16: {
Chris@16: public:
Chris@16:     typedef IntType result_type;
Chris@16:     BOOST_STATIC_CONSTANT(std::size_t, word_size = w);
Chris@16:     BOOST_STATIC_CONSTANT(std::size_t, long_lag = r);
Chris@16:     BOOST_STATIC_CONSTANT(std::size_t, short_lag = s);
Chris@16:     BOOST_STATIC_CONSTANT(uint32_t, default_seed = 19780503u);
Chris@16: 
Chris@16:     // Required by the old Boost.Random concepts
Chris@16:     BOOST_STATIC_CONSTANT(bool, has_fixed_range = false);
Chris@16:     // Backwards compatibility
Chris@16:     BOOST_STATIC_CONSTANT(result_type, modulus = (result_type(1) << w));
Chris@16:     
Chris@16:     BOOST_STATIC_ASSERT(std::numeric_limits<result_type>::is_integer);
Chris@16: 
Chris@16:     /**
Chris@16:      * Constructs a new @c subtract_with_carry_engine and seeds
Chris@16:      * it with the default seed.
Chris@16:      */
Chris@16:     subtract_with_carry_engine() { seed(); }
Chris@16:     /**
Chris@16:      * Constructs a new @c subtract_with_carry_engine and seeds
Chris@16:      * it with @c value.
Chris@16:      */
Chris@16:     BOOST_RANDOM_DETAIL_ARITHMETIC_CONSTRUCTOR(subtract_with_carry_engine,
Chris@16:                                                IntType, value)
Chris@16:     { seed(value); }
Chris@16:     /**
Chris@16:      * Constructs a new @c subtract_with_carry_engine and seeds
Chris@16:      * it with values produced by @c seq.generate().
Chris@16:      */
Chris@16:     BOOST_RANDOM_DETAIL_SEED_SEQ_CONSTRUCTOR(subtract_with_carry_engine,
Chris@16:                                              SeedSeq, seq)
Chris@16:     { seed(seq); }
Chris@16:     /**
Chris@16:      * Constructs a new @c subtract_with_carry_engine and seeds
Chris@16:      * it with values from a range.  first is updated to point
Chris@16:      * one past the last value consumed.  If there are not
Chris@16:      * enough elements in the range to fill the entire state of
Chris@16:      * the generator, throws @c std::invalid_argument.
Chris@16:      */
Chris@16:     template<class It> subtract_with_carry_engine(It& first, It last)
Chris@16:     { seed(first,last); }
Chris@16: 
Chris@16:     // compiler-generated copy ctor and assignment operator are fine
Chris@16: 
Chris@16:     /** Seeds the generator with the default seed. */
Chris@16:     void seed() { seed(default_seed); }
Chris@16:     BOOST_RANDOM_DETAIL_ARITHMETIC_SEED(subtract_with_carry_engine,
Chris@16:                                         IntType, value)
Chris@16:     {
Chris@16:         typedef linear_congruential_engine<uint32_t,40014,0,2147483563> gen_t;
Chris@16:         gen_t intgen(static_cast<boost::uint32_t>(value == 0 ? default_seed : value));
Chris@16:         detail::generator_seed_seq<gen_t> gen(intgen);
Chris@16:         seed(gen);
Chris@16:     }
Chris@16: 
Chris@16:     /** Seeds the generator with values produced by @c seq.generate(). */
Chris@16:     BOOST_RANDOM_DETAIL_SEED_SEQ_SEED(subtract_with_carry, SeedSeq, seq)
Chris@16:     {
Chris@16:         detail::seed_array_int<w>(seq, x);
Chris@16:         carry = (x[long_lag-1] == 0);
Chris@16:         k = 0;
Chris@16:     }
Chris@16: 
Chris@16:     /**
Chris@16:      * Seeds the generator with values from a range.  Updates @c first to
Chris@16:      * point one past the last consumed value.  If the range does not
Chris@16:      * contain enough elements to fill the entire state of the generator,
Chris@16:      * throws @c std::invalid_argument.
Chris@16:      */
Chris@16:     template<class It>
Chris@16:     void seed(It& first, It last)
Chris@16:     {
Chris@16:         detail::fill_array_int<w>(first, last, x);
Chris@16:         carry = (x[long_lag-1] == 0);
Chris@16:         k = 0;
Chris@16:     }
Chris@16: 
Chris@16:     /** Returns the smallest value that the generator can produce. */
Chris@16:     static result_type min BOOST_PREVENT_MACRO_SUBSTITUTION ()
Chris@16:     { return 0; }
Chris@16:     /** Returns the largest value that the generator can produce. */
Chris@16:     static result_type max BOOST_PREVENT_MACRO_SUBSTITUTION ()
Chris@16:     { return boost::low_bits_mask_t<w>::sig_bits; }
Chris@16: 
Chris@16:     /** Returns the next value of the generator. */
Chris@16:     result_type operator()()
Chris@16:     {
Chris@16:         std::size_t short_index =
Chris@16:             (k < short_lag)?
Chris@16:                 (k + long_lag - short_lag) :
Chris@16:                 (k - short_lag);
Chris@16:         carry = do_update(k, short_index, carry);
Chris@16:         IntType result = x[k];
Chris@16:         ++k;
Chris@16:         if(k >= long_lag)
Chris@16:             k = 0;
Chris@16:         return result;
Chris@16:     }
Chris@16: 
Chris@16:     /** Advances the state of the generator by @c z. */
Chris@16:     void discard(boost::uintmax_t z)
Chris@16:     {
Chris@16:         detail::subtract_with_carry_discard::apply(*this, z);
Chris@16:     }
Chris@16: 
Chris@16:     /** Fills a range with random values. */
Chris@16:     template<class It>
Chris@16:     void generate(It first, It last)
Chris@16:     { detail::generate_from_int(*this, first, last); }
Chris@16:  
Chris@16:     /** Writes a @c subtract_with_carry_engine to a @c std::ostream. */
Chris@16:     BOOST_RANDOM_DETAIL_OSTREAM_OPERATOR(os, subtract_with_carry_engine, f)
Chris@16:     {
Chris@16:         for(unsigned int j = 0; j < f.long_lag; ++j)
Chris@16:             os << f.compute(j) << ' ';
Chris@16:         os << f.carry;
Chris@16:         return os;
Chris@16:     }
Chris@16: 
Chris@16:     /** Reads a @c subtract_with_carry_engine from a @c std::istream. */
Chris@16:     BOOST_RANDOM_DETAIL_ISTREAM_OPERATOR(is, subtract_with_carry_engine, f)
Chris@16:     {
Chris@16:         for(unsigned int j = 0; j < f.long_lag; ++j)
Chris@16:             is >> f.x[j] >> std::ws;
Chris@16:         is >> f.carry;
Chris@16:         f.k = 0;
Chris@16:         return is;
Chris@16:     }
Chris@16: 
Chris@16:     /**
Chris@16:      * Returns true if the two generators will produce identical
Chris@16:      * sequences of values.
Chris@16:      */
Chris@16:     BOOST_RANDOM_DETAIL_EQUALITY_OPERATOR(subtract_with_carry_engine, x, y)
Chris@16:     {
Chris@16:         for(unsigned int j = 0; j < r; ++j)
Chris@16:             if(x.compute(j) != y.compute(j))
Chris@16:                 return false;
Chris@16:         return true;
Chris@16:     }
Chris@16: 
Chris@16:     /**
Chris@16:      * Returns true if the two generators will produce different
Chris@16:      * sequences of values.
Chris@16:      */
Chris@16:     BOOST_RANDOM_DETAIL_INEQUALITY_OPERATOR(subtract_with_carry_engine)
Chris@16: 
Chris@16: private:
Chris@16:     /// \cond show_private
Chris@16:     // returns x(i-r+index), where index is in 0..r-1
Chris@16:     IntType compute(unsigned int index) const
Chris@16:     {
Chris@16:         return x[(k+index) % long_lag];
Chris@16:     }
Chris@16: 
Chris@16:     friend struct detail::subtract_with_carry_discard;
Chris@16: 
Chris@16:     IntType do_update(std::size_t current, std::size_t short_index, IntType carry)
Chris@16:     {
Chris@16:         IntType delta;
Chris@16:         IntType temp = x[current] + carry;
Chris@16:         if (x[short_index] >= temp) {
Chris@16:             // x(n) >= 0
Chris@16:             delta =  x[short_index] - temp;
Chris@16:             carry = 0;
Chris@16:         } else {
Chris@16:             // x(n) < 0
Chris@16:             delta = modulus - temp + x[short_index];
Chris@16:             carry = 1;
Chris@16:         }
Chris@16:         x[current] = delta;
Chris@16:         return carry;
Chris@16:     }
Chris@16:     /// \endcond
Chris@16: 
Chris@16:     // state representation; next output (state) is x(i)
Chris@16:     //   x[0]  ... x[k] x[k+1] ... x[long_lag-1]     represents
Chris@16:     //  x(i-k) ... x(i) x(i+1) ... x(i-k+long_lag-1)
Chris@16:     // speed: base: 20-25 nsec
Chris@16:     // ranlux_4: 230 nsec, ranlux_7: 430 nsec, ranlux_14: 810 nsec
Chris@16:     // This state representation makes operator== and save/restore more
Chris@16:     // difficult, because we've already computed "too much" and thus
Chris@16:     // have to undo some steps to get at x(i-r) etc.
Chris@16: 
Chris@16:     // state representation: next output (state) is x(i)
Chris@16:     //   x[0]  ... x[k] x[k+1]          ... x[long_lag-1]     represents
Chris@16:     //  x(i-k) ... x(i) x(i-long_lag+1) ... x(i-k-1)
Chris@16:     // speed: base 28 nsec
Chris@16:     // ranlux_4: 370 nsec, ranlux_7: 688 nsec, ranlux_14: 1343 nsec
Chris@16:     IntType x[long_lag];
Chris@16:     std::size_t k;
Chris@16:     IntType carry;
Chris@16: };
Chris@16: 
Chris@16: #ifndef BOOST_NO_INCLASS_MEMBER_INITIALIZATION
Chris@16: //  A definition is required even for integral static constants
Chris@16: template<class IntType, std::size_t w, std::size_t s, std::size_t r>
Chris@16: const bool subtract_with_carry_engine<IntType, w, s, r>::has_fixed_range;
Chris@16: template<class IntType, std::size_t w, std::size_t s, std::size_t r>
Chris@16: const IntType subtract_with_carry_engine<IntType, w, s, r>::modulus;
Chris@16: template<class IntType, std::size_t w, std::size_t s, std::size_t r>
Chris@16: const std::size_t subtract_with_carry_engine<IntType, w, s, r>::word_size;
Chris@16: template<class IntType, std::size_t w, std::size_t s, std::size_t r>
Chris@16: const std::size_t subtract_with_carry_engine<IntType, w, s, r>::long_lag;
Chris@16: template<class IntType, std::size_t w, std::size_t s, std::size_t r>
Chris@16: const std::size_t subtract_with_carry_engine<IntType, w, s, r>::short_lag;
Chris@16: template<class IntType, std::size_t w, std::size_t s, std::size_t r>
Chris@16: const uint32_t subtract_with_carry_engine<IntType, w, s, r>::default_seed;
Chris@16: #endif
Chris@16: 
Chris@16: 
Chris@16: // use a floating-point representation to produce values in [0..1)
Chris@16: /**
Chris@16:  * Instantiations of \subtract_with_carry_01_engine model a
Chris@16:  * \pseudo_random_number_generator.  The algorithm is
Chris@16:  * described in
Chris@16:  *
Chris@16:  *  @blockquote
Chris@16:  *  "A New Class of Random Number Generators", George
Chris@16:  *  Marsaglia and Arif Zaman, Annals of Applied Probability,
Chris@16:  *  Volume 1, Number 3 (1991), 462-480.
Chris@16:  *  @endblockquote
Chris@16:  */
Chris@16: template<class RealType, std::size_t w, std::size_t s, std::size_t r>
Chris@16: class subtract_with_carry_01_engine
Chris@16: {
Chris@16: public:
Chris@16:     typedef RealType result_type;
Chris@16:     BOOST_STATIC_CONSTANT(bool, has_fixed_range = false);
Chris@16:     BOOST_STATIC_CONSTANT(std::size_t, word_size = w);
Chris@16:     BOOST_STATIC_CONSTANT(std::size_t, long_lag = r);
Chris@16:     BOOST_STATIC_CONSTANT(std::size_t, short_lag = s);
Chris@16:     BOOST_STATIC_CONSTANT(boost::uint32_t, default_seed = 19780503u);
Chris@16: 
Chris@16:     BOOST_STATIC_ASSERT(!std::numeric_limits<result_type>::is_integer);
Chris@16: 
Chris@16:     /** Creates a new \subtract_with_carry_01_engine using the default seed. */
Chris@16:     subtract_with_carry_01_engine() { init_modulus(); seed(); }
Chris@16:     /** Creates a new subtract_with_carry_01_engine and seeds it with value. */
Chris@16:     BOOST_RANDOM_DETAIL_ARITHMETIC_CONSTRUCTOR(subtract_with_carry_01_engine,
Chris@16:                                                boost::uint32_t, value)
Chris@16:     { init_modulus(); seed(value); }
Chris@16:     /**
Chris@16:      * Creates a new \subtract_with_carry_01_engine and seeds with values
Chris@16:      * produced by seq.generate().
Chris@16:      */
Chris@16:     BOOST_RANDOM_DETAIL_SEED_SEQ_CONSTRUCTOR(subtract_with_carry_01_engine,
Chris@16:                                              SeedSeq, seq)
Chris@16:     { init_modulus(); seed(seq); }
Chris@16:     /**
Chris@16:      * Creates a new \subtract_with_carry_01_engine and seeds it with values
Chris@16:      * from a range.  Advances first to point one past the last consumed
Chris@16:      * value.  If the range does not contain enough elements to fill the
Chris@16:      * entire state, throws @c std::invalid_argument.
Chris@16:      */
Chris@16:     template<class It> subtract_with_carry_01_engine(It& first, It last)
Chris@16:     { init_modulus(); seed(first,last); }
Chris@16: 
Chris@16: private:
Chris@16:     /// \cond show_private
Chris@16:     void init_modulus()
Chris@16:     {
Chris@16: #ifndef BOOST_NO_STDC_NAMESPACE
Chris@16:         // allow for Koenig lookup
Chris@16:         using std::pow;
Chris@16: #endif
Chris@16:         _modulus = pow(RealType(2), RealType(word_size));
Chris@16:     }
Chris@16:     /// \endcond
Chris@16: 
Chris@16: public:
Chris@16:     // compiler-generated copy ctor and assignment operator are fine
Chris@16: 
Chris@16:     /** Seeds the generator with the default seed. */
Chris@16:     void seed() { seed(default_seed); }
Chris@16: 
Chris@16:     /** Seeds the generator with @c value. */
Chris@16:     BOOST_RANDOM_DETAIL_ARITHMETIC_SEED(subtract_with_carry_01_engine,
Chris@16:                                         boost::uint32_t, value)
Chris@16:     {
Chris@16:         typedef linear_congruential_engine<uint32_t, 40014, 0, 2147483563> gen_t;
Chris@16:         gen_t intgen(value == 0 ? default_seed : value);
Chris@16:         detail::generator_seed_seq<gen_t> gen(intgen);
Chris@16:         seed(gen);
Chris@16:     }
Chris@16: 
Chris@16:     /** Seeds the generator with values produced by @c seq.generate(). */
Chris@16:     BOOST_RANDOM_DETAIL_SEED_SEQ_SEED(subtract_with_carry_01_engine,
Chris@16:                                       SeedSeq, seq)
Chris@16:     {
Chris@16:         detail::seed_array_real<w>(seq, x);
Chris@16:         carry = (x[long_lag-1] ? 0 : 1 / _modulus);
Chris@16:         k = 0;
Chris@16:     }
Chris@16: 
Chris@16:     /**
Chris@16:      * Seeds the generator with values from a range.  Updates first to
Chris@16:      * point one past the last consumed element.  If there are not
Chris@16:      * enough elements in the range to fill the entire state, throws
Chris@16:      * @c std::invalid_argument.
Chris@16:      */
Chris@16:     template<class It>
Chris@16:     void seed(It& first, It last)
Chris@16:     {
Chris@16:         detail::fill_array_real<w>(first, last, x);
Chris@16:         carry = (x[long_lag-1] ? 0 : 1 / _modulus);
Chris@16:         k = 0;
Chris@16:     }
Chris@16: 
Chris@16:     /** Returns the smallest value that the generator can produce. */
Chris@16:     static result_type min BOOST_PREVENT_MACRO_SUBSTITUTION ()
Chris@16:     { return result_type(0); }
Chris@16:     /** Returns the largest value that the generator can produce. */
Chris@16:     static result_type max BOOST_PREVENT_MACRO_SUBSTITUTION ()
Chris@16:     { return result_type(1); }
Chris@16: 
Chris@16:     /** Returns the next value of the generator. */
Chris@16:     result_type operator()()
Chris@16:     {
Chris@16:         std::size_t short_index =
Chris@16:             (k < short_lag) ?
Chris@16:                 (k + long_lag - short_lag) :
Chris@16:                 (k - short_lag);
Chris@16:         carry = do_update(k, short_index, carry);
Chris@16:         RealType result = x[k];
Chris@16:         ++k;
Chris@16:         if(k >= long_lag)
Chris@16:             k = 0;
Chris@16:         return result;
Chris@16:     }
Chris@16: 
Chris@16:     /** Advances the state of the generator by @c z. */
Chris@16:     void discard(boost::uintmax_t z)
Chris@16:     { detail::subtract_with_carry_discard::apply(*this, z); }
Chris@16: 
Chris@16:     /** Fills a range with random values. */
Chris@16:     template<class Iter>
Chris@16:     void generate(Iter first, Iter last)
Chris@16:     { detail::generate_from_real(*this, first, last); }
Chris@16: 
Chris@16:     /** Writes a \subtract_with_carry_01_engine to a @c std::ostream. */
Chris@16:     BOOST_RANDOM_DETAIL_OSTREAM_OPERATOR(os, subtract_with_carry_01_engine, f)
Chris@16:     {
Chris@16:         std::ios_base::fmtflags oldflags =
Chris@16:             os.flags(os.dec | os.fixed | os.left); 
Chris@16:         for(unsigned int j = 0; j < f.long_lag; ++j)
Chris@16:             os << (f.compute(j) * f._modulus) << ' ';
Chris@16:         os << (f.carry * f._modulus);
Chris@16:         os.flags(oldflags);
Chris@16:         return os;
Chris@16:     }
Chris@16:     
Chris@16:     /** Reads a \subtract_with_carry_01_engine from a @c std::istream. */
Chris@16:     BOOST_RANDOM_DETAIL_ISTREAM_OPERATOR(is, subtract_with_carry_01_engine, f)
Chris@16:     {
Chris@16:         RealType value;
Chris@16:         for(unsigned int j = 0; j < long_lag; ++j) {
Chris@16:             is >> value >> std::ws;
Chris@16:             f.x[j] = value / f._modulus;
Chris@16:         }
Chris@16:         is >> value;
Chris@16:         f.carry = value / f._modulus;
Chris@16:         f.k = 0;
Chris@16:         return is;
Chris@16:     }
Chris@16: 
Chris@16:     /** Returns true if the two generators will produce identical sequences. */
Chris@16:     BOOST_RANDOM_DETAIL_EQUALITY_OPERATOR(subtract_with_carry_01_engine, x, y)
Chris@16:     {
Chris@16:         for(unsigned int j = 0; j < r; ++j)
Chris@16:             if(x.compute(j) != y.compute(j))
Chris@16:                 return false;
Chris@16:         return true;
Chris@16:     }
Chris@16: 
Chris@16:     /** Returns true if the two generators will produce different sequences. */
Chris@16:     BOOST_RANDOM_DETAIL_INEQUALITY_OPERATOR(subtract_with_carry_01_engine)
Chris@16: 
Chris@16: private:
Chris@16:     /// \cond show_private
Chris@16:     RealType compute(unsigned int index) const
Chris@16:     {
Chris@16:         return x[(k+index) % long_lag];
Chris@16:     }
Chris@16: 
Chris@16:     friend struct detail::subtract_with_carry_discard;
Chris@16: 
Chris@16:     RealType do_update(std::size_t current, std::size_t short_index, RealType carry)
Chris@16:     {
Chris@16:         RealType delta = x[short_index] - x[current] - carry;
Chris@16:         if(delta < 0) {
Chris@16:             delta += RealType(1);
Chris@16:             carry = RealType(1)/_modulus;
Chris@16:         } else {
Chris@16:             carry = 0;
Chris@16:         }
Chris@16:         x[current] = delta;
Chris@16:         return carry;
Chris@16:     }
Chris@16:     /// \endcond
Chris@16:     std::size_t k;
Chris@16:     RealType carry;
Chris@16:     RealType x[long_lag];
Chris@16:     RealType _modulus;
Chris@16: };
Chris@16: 
Chris@16: #ifndef BOOST_NO_INCLASS_MEMBER_INITIALIZATION
Chris@16: //  A definition is required even for integral static constants
Chris@16: template<class RealType, std::size_t w, std::size_t s, std::size_t r>
Chris@16: const bool subtract_with_carry_01_engine<RealType, w, s, r>::has_fixed_range;
Chris@16: template<class RealType, std::size_t w, std::size_t s, std::size_t r>
Chris@16: const std::size_t subtract_with_carry_01_engine<RealType, w, s, r>::word_size;
Chris@16: template<class RealType, std::size_t w, std::size_t s, std::size_t r>
Chris@16: const std::size_t subtract_with_carry_01_engine<RealType, w, s, r>::long_lag;
Chris@16: template<class RealType, std::size_t w, std::size_t s, std::size_t r>
Chris@16: const std::size_t subtract_with_carry_01_engine<RealType, w, s, r>::short_lag;
Chris@16: template<class RealType, std::size_t w, std::size_t s, std::size_t r>
Chris@16: const uint32_t subtract_with_carry_01_engine<RealType, w, s, r>::default_seed;
Chris@16: #endif
Chris@16: 
Chris@16: 
Chris@16: /// \cond show_deprecated
Chris@16: 
Chris@16: template<class IntType, IntType m, unsigned s, unsigned r, IntType v>
Chris@16: class subtract_with_carry :
Chris@16:     public subtract_with_carry_engine<IntType,
Chris@16:         boost::static_log2<m>::value, s, r>
Chris@16: {
Chris@16:     typedef subtract_with_carry_engine<IntType,
Chris@16:         boost::static_log2<m>::value, s, r> base_type;
Chris@16: public:
Chris@16:     subtract_with_carry() {}
Chris@16:     BOOST_RANDOM_DETAIL_GENERATOR_CONSTRUCTOR(subtract_with_carry, Gen, gen)
Chris@16:     { seed(gen); }
Chris@16:     BOOST_RANDOM_DETAIL_ARITHMETIC_CONSTRUCTOR(subtract_with_carry,
Chris@16:                                                IntType, val)
Chris@16:     { seed(val); }
Chris@16:     template<class It>
Chris@16:     subtract_with_carry(It& first, It last) : base_type(first, last) {}
Chris@16:     void seed() { base_type::seed(); }
Chris@16:     BOOST_RANDOM_DETAIL_GENERATOR_SEED(subtract_with_carry, Gen, gen)
Chris@16:     {
Chris@16:         detail::generator_seed_seq<Gen> seq(gen);
Chris@16:         base_type::seed(seq);
Chris@16:     }
Chris@16:     BOOST_RANDOM_DETAIL_ARITHMETIC_SEED(subtract_with_carry, IntType, val)
Chris@16:     { base_type::seed(val); }
Chris@16:     template<class It>
Chris@16:     void seed(It& first, It last) { base_type::seed(first, last); }
Chris@16: };
Chris@16: 
Chris@16: template<class RealType, int w, unsigned s, unsigned r, int v = 0>
Chris@16: class subtract_with_carry_01 :
Chris@16:     public subtract_with_carry_01_engine<RealType, w, s, r>
Chris@16: {
Chris@16:     typedef subtract_with_carry_01_engine<RealType, w, s, r> base_type;
Chris@16: public:
Chris@16:     subtract_with_carry_01() {}
Chris@16:     BOOST_RANDOM_DETAIL_GENERATOR_CONSTRUCTOR(subtract_with_carry_01, Gen, gen)
Chris@16:     { seed(gen); }
Chris@16:     BOOST_RANDOM_DETAIL_ARITHMETIC_CONSTRUCTOR(subtract_with_carry_01,
Chris@16:                                                uint32_t, val)
Chris@16:     { seed(val); }
Chris@16:     template<class It>
Chris@16:     subtract_with_carry_01(It& first, It last) : base_type(first, last) {}
Chris@16:     void seed() { base_type::seed(); }
Chris@16:     BOOST_RANDOM_DETAIL_GENERATOR_SEED(subtract_with_carry_01, Gen, gen)
Chris@16:     {
Chris@16:         detail::generator_seed_seq<Gen> seq(gen);
Chris@16:         base_type::seed(seq);
Chris@16:     }
Chris@16:     BOOST_RANDOM_DETAIL_ARITHMETIC_SEED(subtract_with_carry_01, uint32_t, val)
Chris@16:     { base_type::seed(val); }
Chris@16:     template<class It>
Chris@16:     void seed(It& first, It last) { base_type::seed(first, last); }
Chris@16: };
Chris@16: 
Chris@16: /// \endcond
Chris@16: 
Chris@16: namespace detail {
Chris@16: 
Chris@16: template<class Engine>
Chris@16: struct generator_bits;
Chris@16: 
Chris@16: template<class RealType, std::size_t w, std::size_t s, std::size_t r>
Chris@16: struct generator_bits<subtract_with_carry_01_engine<RealType, w, s, r> > {
Chris@16:     static std::size_t value() { return w; }
Chris@16: };
Chris@16: 
Chris@16: template<class RealType, int w, unsigned s, unsigned r, int v>
Chris@16: struct generator_bits<subtract_with_carry_01<RealType, w, s, r, v> > {
Chris@16:     static std::size_t value() { return w; }
Chris@16: };
Chris@16: 
Chris@16: }
Chris@16: 
Chris@16: } // namespace random
Chris@16: } // namespace boost
Chris@16: 
Chris@16: #endif // BOOST_RANDOM_SUBTRACT_WITH_CARRY_HPP