1 /* boost random/inversive_congruential.hpp header file 2 * 3 * Copyright Jens Maurer 2000-2001 4 * Distributed under the Boost Software License, Version 1.0. (See 5 * accompanying file LICENSE_1_0.txt or copy at 6 * http://www.boost.org/LICENSE_1_0.txt) 7 * 8 * See http://www.boost.org for most recent version including documentation. 9 * 10 * $Id$ 11 * 12 * Revision history 13 * 2001-02-18 moved to individual header files 14 */ 15 16 #ifndef BOOST_RANDOM_INVERSIVE_CONGRUENTIAL_HPP 17 #define BOOST_RANDOM_INVERSIVE_CONGRUENTIAL_HPP 18 19 #include <iosfwd> 20 #include <stdexcept> 21 #include <boost/assert.hpp> 22 #include <boost/config.hpp> 23 #include <boost/cstdint.hpp> 24 #include <boost/integer/static_log2.hpp> 25 #include <boost/random/detail/config.hpp> 26 #include <boost/random/detail/const_mod.hpp> 27 #include <boost/random/detail/seed.hpp> 28 #include <boost/random/detail/operators.hpp> 29 #include <boost/random/detail/seed_impl.hpp> 30 31 #include <boost/random/detail/disable_warnings.hpp> 32 33 namespace boost { 34 namespace random { 35 36 // Eichenauer and Lehn 1986 37 /** 38 * Instantiations of class template @c inversive_congruential_engine model a 39 * \pseudo_random_number_generator. It uses the inversive congruential 40 * algorithm (ICG) described in 41 * 42 * @blockquote 43 * "Inversive pseudorandom number generators: concepts, results and links", 44 * Peter Hellekalek, In: "Proceedings of the 1995 Winter Simulation 45 * Conference", C. Alexopoulos, K. Kang, W.R. Lilegdon, and D. Goldsman 46 * (editors), 1995, pp. 255-262. ftp://random.mat.sbg.ac.at/pub/data/wsc95.ps 47 * @endblockquote 48 * 49 * The output sequence is defined by x(n+1) = (a*inv(x(n)) - b) (mod p), 50 * where x(0), a, b, and the prime number p are parameters of the generator. 51 * The expression inv(k) denotes the multiplicative inverse of k in the 52 * field of integer numbers modulo p, with inv(0) := 0. 53 * 54 * The template parameter IntType shall denote a signed integral type large 55 * enough to hold p; a, b, and p are the parameters of the generators. The 56 * template parameter val is the validation value checked by validation. 57 * 58 * @xmlnote 59 * The implementation currently uses the Euclidian Algorithm to compute 60 * the multiplicative inverse. Therefore, the inversive generators are about 61 * 10-20 times slower than the others (see section"performance"). However, 62 * the paper talks of only 3x slowdown, so the Euclidian Algorithm is probably 63 * not optimal for calculating the multiplicative inverse. 64 * @endxmlnote 65 */ 66 template<class IntType, IntType a, IntType b, IntType p> 67 class inversive_congruential_engine 68 { 69 public: 70 typedef IntType result_type; 71 BOOST_STATIC_CONSTANT(bool, has_fixed_range = false); 72 73 BOOST_STATIC_CONSTANT(result_type, multiplier = a); 74 BOOST_STATIC_CONSTANT(result_type, increment = b); 75 BOOST_STATIC_CONSTANT(result_type, modulus = p); 76 BOOST_STATIC_CONSTANT(IntType, default_seed = 1); 77 BOOST_PREVENT_MACRO_SUBSTITUTION()78 static result_type min BOOST_PREVENT_MACRO_SUBSTITUTION () { return b == 0 ? 1 : 0; } BOOST_PREVENT_MACRO_SUBSTITUTION()79 static result_type max BOOST_PREVENT_MACRO_SUBSTITUTION () { return p-1; } 80 81 /** 82 * Constructs an @c inversive_congruential_engine, seeding it with 83 * the default seed. 84 */ inversive_congruential_engine()85 inversive_congruential_engine() { seed(); } 86 87 /** 88 * Constructs an @c inversive_congruential_engine, seeding it with @c x0. 89 */ BOOST_RANDOM_DETAIL_ARITHMETIC_CONSTRUCTOR(inversive_congruential_engine,IntType,x0)90 BOOST_RANDOM_DETAIL_ARITHMETIC_CONSTRUCTOR(inversive_congruential_engine, 91 IntType, x0) 92 { seed(x0); } 93 94 /** 95 * Constructs an @c inversive_congruential_engine, seeding it with values 96 * produced by a call to @c seq.generate(). 97 */ BOOST_RANDOM_DETAIL_SEED_SEQ_CONSTRUCTOR(inversive_congruential_engine,SeedSeq,seq)98 BOOST_RANDOM_DETAIL_SEED_SEQ_CONSTRUCTOR(inversive_congruential_engine, 99 SeedSeq, seq) 100 { seed(seq); } 101 102 /** 103 * Constructs an @c inversive_congruential_engine, seeds it 104 * with values taken from the itrator range [first, last), 105 * and adjusts first to point to the element after the last one 106 * used. If there are not enough elements, throws @c std::invalid_argument. 107 * 108 * first and last must be input iterators. 109 */ inversive_congruential_engine(It & first,It last)110 template<class It> inversive_congruential_engine(It& first, It last) 111 { seed(first, last); } 112 113 /** 114 * Calls seed(default_seed) 115 */ seed()116 void seed() { seed(default_seed); } 117 118 /** 119 * If c mod m is zero and x0 mod m is zero, changes the current value of 120 * the generator to 1. Otherwise, changes it to x0 mod m. If c is zero, 121 * distinct seeds in the range [1,m) will leave the generator in distinct 122 * states. If c is not zero, the range is [0,m). 123 */ BOOST_RANDOM_DETAIL_ARITHMETIC_SEED(inversive_congruential_engine,IntType,x0)124 BOOST_RANDOM_DETAIL_ARITHMETIC_SEED(inversive_congruential_engine, IntType, x0) 125 { 126 // wrap _x if it doesn't fit in the destination 127 if(modulus == 0) { 128 _value = x0; 129 } else { 130 _value = x0 % modulus; 131 } 132 // handle negative seeds 133 if(_value <= 0 && _value != 0) { 134 _value += modulus; 135 } 136 // adjust to the correct range 137 if(increment == 0 && _value == 0) { 138 _value = 1; 139 } 140 BOOST_ASSERT(_value >= (min)()); 141 BOOST_ASSERT(_value <= (max)()); 142 } 143 144 /** 145 * Seeds an @c inversive_congruential_engine using values from a SeedSeq. 146 */ BOOST_RANDOM_DETAIL_SEED_SEQ_SEED(inversive_congruential_engine,SeedSeq,seq)147 BOOST_RANDOM_DETAIL_SEED_SEQ_SEED(inversive_congruential_engine, SeedSeq, seq) 148 { seed(detail::seed_one_int<IntType, modulus>(seq)); } 149 150 /** 151 * seeds an @c inversive_congruential_engine with values taken 152 * from the itrator range [first, last) and adjusts @c first to 153 * point to the element after the last one used. If there are 154 * not enough elements, throws @c std::invalid_argument. 155 * 156 * @c first and @c last must be input iterators. 157 */ seed(It & first,It last)158 template<class It> void seed(It& first, It last) 159 { seed(detail::get_one_int<IntType, modulus>(first, last)); } 160 161 /** Returns the next output of the generator. */ operator ()()162 IntType operator()() 163 { 164 typedef const_mod<IntType, p> do_mod; 165 _value = do_mod::mult_add(a, do_mod::invert(_value), b); 166 return _value; 167 } 168 169 /** Fills a range with random values */ 170 template<class Iter> generate(Iter first,Iter last)171 void generate(Iter first, Iter last) 172 { detail::generate_from_int(*this, first, last); } 173 174 /** Advances the state of the generator by @c z. */ discard(boost::uintmax_t z)175 void discard(boost::uintmax_t z) 176 { 177 for(boost::uintmax_t j = 0; j < z; ++j) { 178 (*this)(); 179 } 180 } 181 182 /** 183 * Writes the textual representation of the generator to a @c std::ostream. 184 */ BOOST_RANDOM_DETAIL_OSTREAM_OPERATOR(os,inversive_congruential_engine,x)185 BOOST_RANDOM_DETAIL_OSTREAM_OPERATOR(os, inversive_congruential_engine, x) 186 { 187 os << x._value; 188 return os; 189 } 190 191 /** 192 * Reads the textual representation of the generator from a @c std::istream. 193 */ BOOST_RANDOM_DETAIL_ISTREAM_OPERATOR(is,inversive_congruential_engine,x)194 BOOST_RANDOM_DETAIL_ISTREAM_OPERATOR(is, inversive_congruential_engine, x) 195 { 196 is >> x._value; 197 return is; 198 } 199 200 /** 201 * Returns true if the two generators will produce identical 202 * sequences of outputs. 203 */ BOOST_RANDOM_DETAIL_EQUALITY_OPERATOR(inversive_congruential_engine,x,y)204 BOOST_RANDOM_DETAIL_EQUALITY_OPERATOR(inversive_congruential_engine, x, y) 205 { return x._value == y._value; } 206 207 /** 208 * Returns true if the two generators will produce different 209 * sequences of outputs. 210 */ 211 BOOST_RANDOM_DETAIL_INEQUALITY_OPERATOR(inversive_congruential_engine) 212 213 private: 214 IntType _value; 215 }; 216 217 #ifndef BOOST_NO_INCLASS_MEMBER_INITIALIZATION 218 // A definition is required even for integral static constants 219 template<class IntType, IntType a, IntType b, IntType p> 220 const bool inversive_congruential_engine<IntType, a, b, p>::has_fixed_range; 221 template<class IntType, IntType a, IntType b, IntType p> 222 const typename inversive_congruential_engine<IntType, a, b, p>::result_type inversive_congruential_engine<IntType, a, b, p>::multiplier; 223 template<class IntType, IntType a, IntType b, IntType p> 224 const typename inversive_congruential_engine<IntType, a, b, p>::result_type inversive_congruential_engine<IntType, a, b, p>::increment; 225 template<class IntType, IntType a, IntType b, IntType p> 226 const typename inversive_congruential_engine<IntType, a, b, p>::result_type inversive_congruential_engine<IntType, a, b, p>::modulus; 227 template<class IntType, IntType a, IntType b, IntType p> 228 const typename inversive_congruential_engine<IntType, a, b, p>::result_type inversive_congruential_engine<IntType, a, b, p>::default_seed; 229 #endif 230 231 /// \cond show_deprecated 232 233 // provided for backwards compatibility 234 template<class IntType, IntType a, IntType b, IntType p, IntType val = 0> 235 class inversive_congruential : public inversive_congruential_engine<IntType, a, b, p> 236 { 237 typedef inversive_congruential_engine<IntType, a, b, p> base_type; 238 public: inversive_congruential(IntType x0=1)239 inversive_congruential(IntType x0 = 1) : base_type(x0) {} 240 template<class It> inversive_congruential(It & first,It last)241 inversive_congruential(It& first, It last) : base_type(first, last) {} 242 }; 243 244 /// \endcond 245 246 /** 247 * The specialization hellekalek1995 was suggested in 248 * 249 * @blockquote 250 * "Inversive pseudorandom number generators: concepts, results and links", 251 * Peter Hellekalek, In: "Proceedings of the 1995 Winter Simulation 252 * Conference", C. Alexopoulos, K. Kang, W.R. Lilegdon, and D. Goldsman 253 * (editors), 1995, pp. 255-262. ftp://random.mat.sbg.ac.at/pub/data/wsc95.ps 254 * @endblockquote 255 */ 256 typedef inversive_congruential_engine<uint32_t, 9102, 2147483647-36884165, 257 2147483647> hellekalek1995; 258 259 } // namespace random 260 261 using random::hellekalek1995; 262 263 } // namespace boost 264 265 #include <boost/random/detail/enable_warnings.hpp> 266 267 #endif // BOOST_RANDOM_INVERSIVE_CONGRUENTIAL_HPP 268