Horizon
random.hpp
Go to the documentation of this file.
1 // Range v3 library
3 //
4 // Copyright Casey Carter 2016
5 //
6 // Use, modification and distribution is subject to the
7 // Boost Software License, Version 1.0. (See accompanying
8 // file LICENSE_1_0.txt or copy at
9 // http://www.boost.org/LICENSE_1_0.txt)
10 //
11 // Project home: https://github.com/ericniebler/range-v3
12 //
13 
14 /*
15  * Random-Number Utilities (randutil)
16  * Addresses common issues with C++11 random number generation.
17  * Makes good seeding easier, and makes using RNGs easy while retaining
18  * all the power.
19  *
20  * The MIT License (MIT)
21  *
22  * Copyright (c) 2015 Melissa E. O'Neill
23  *
24  * Permission is hereby granted, free of charge, to any person obtaining a copy
25  * of this software and associated documentation files (the "Software"), to deal
26  * in the Software without restriction, including without limitation the rights
27  * to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
28  * copies of the Software, and to permit persons to whom the Software is
29  * furnished to do so, subject to the following conditions:
30  *
31  * The above copyright notice and this permission notice shall be included in
32  * all copies or substantial portions of the Software.
33  *
34  * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
35  * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
36  * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
37  * AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
38  * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
39  * OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
40  * SOFTWARE.
41  */
42 
43 #ifndef RANGES_V3_UTILITY_RANDOM_HPP
44 #define RANGES_V3_UTILITY_RANDOM_HPP
45 
46 #include <array>
47 #include <cstddef>
48 #include <cstdint>
49 #include <initializer_list>
50 #include <new>
51 #include <random>
52 
53 #include <meta/meta.hpp>
54 
55 #include <concepts/concepts.hpp>
56 
57 #include <range/v3/range_fwd.hpp>
58 
64 
65 #if !RANGES_CXX_THREAD_LOCAL
66 #include <mutex>
67 #endif
68 
69 #include <range/v3/detail/prologue.hpp>
70 
71 RANGES_DIAGNOSTIC_PUSH
72 RANGES_DIAGNOSTIC_IGNORE_CXX17_COMPAT
73 
74 namespace ranges
75 {
78  // clang-format off
81  template<typename Gen>
82  CPP_requires(uniform_random_bit_generator_,
83  requires() //
84  (
85  Gen::min(),
86  Gen::max()
87  ));
90  template(typename Gen)(
91  concept (uniform_random_bit_generator_)(Gen),
92  unsigned_integral<invoke_result_t<Gen &>> AND
93  same_as<invoke_result_t<Gen &>, decltype(Gen::min())> AND
94  same_as<invoke_result_t<Gen &>, decltype(Gen::max())>);
95 
98  template<typename Gen>
100  invocable<Gen &> &&
101  CPP_requires_ref(ranges::uniform_random_bit_generator_, Gen) &&
102  CPP_concept_ref(ranges::uniform_random_bit_generator_, Gen);
103  // clang-format on
105 
107  namespace detail
108  {
109  namespace randutils
110  {
111  inline std::array<std::uint32_t, 8> get_entropy()
112  {
113  std::array<std::uint32_t, 8> seeds;
114 
115  // Hopefully high-quality entropy from random_device.
116 #if defined(__GLIBCXX__) && defined(RANGES_WORKAROUND_VALGRIND_RDRAND)
117  std::random_device rd{"/dev/urandom"};
118 #else
119  std::random_device rd;
120 #endif
121  std::uniform_int_distribution<std::uint32_t> dist{};
122  ranges::generate(seeds, [&] { return dist(rd); });
123 
124  return seeds;
125  }
126 
127  template(typename I)(
128  requires unsigned_integral<I>)
129  constexpr I fast_exp(I x, I power, I result = I{1})
130  {
131  return power == I{0}
132  ? result
133  : randutils::fast_exp(
134  x * x, power >> 1, result * (power & I{1} ? x : 1));
135  }
136 
138  //
139  // seed_seq_fe
140  //
142 
143  /*
144  * seed_seq_fe implements a fixed-entropy seed sequence; it conforms to all
145  * the requirements of a Seed Sequence concept.
146  *
147  * seed_seq_fe<N> implements a seed sequence which seeds based on a store of
148  * N * 32 bits of entropy. Typically, it would be initialized with N or more
149  * integers.
150  *
151  * seed_seq_fe128 and seed_seq_fe256 are provided as convenience typedefs for
152  * 128- and 256-bit entropy stores respectively. These variants outperform
153  * std::seed_seq, while being better mixing the bits it is provided as
154  * entropy. In almost all common use cases, they serve as better drop-in
155  * replacements for seed_seq.
156  *
157  * Technical details
158  *
159  * Assuming it constructed with M seed integers as input, it exhibits the
160  * following properties
161  *
162  * * Diffusion/Avalanche: A single-bit change in any of the M inputs has a
163  * 50% chance of flipping every bit in the bitstream produced by generate.
164  * Initializing the N-word entropy store with M words requires O(N * M)
165  * time precisely because of the avalanche requirements. Once constructed,
166  * calls to generate are linear in the number of words generated.
167  *
168  * * Bias freedom/Bijection: If M == N, the state of the entropy store is a
169  * bijection from the M inputs (i.e., no states occur twice, none are
170  * omitted). If M > N the number of times each state can occur is the same
171  * (each state occurs 2**(32*(M-N)) times, where ** is the power function).
172  * If M < N, some states cannot occur (bias) but no state occurs more
173  * than once (it's impossible to avoid bias if M < N; ideally N should not
174  * be chosen so that it is more than M).
175  *
176  * Likewise, the generate function has similar properties (with the entropy
177  * store as the input data). If more outputs are requested than there is
178  * entropy, some outputs cannot occur. For example, the Mersenne Twister
179  * will request 624 outputs, to initialize its 19937-bit state, which is
180  * much larger than a 128-bit or 256-bit entropy pool. But in practice,
181  * limiting the Mersenne Twister to 2**128 possible initializations gives
182  * us enough initializations to give a unique initialization to trillions
183  * of computers for billions of years. If you really have 624 words of
184  * *real* high-quality entropy you want to use, you probably don't need
185  * an entropy mixer like this class at all. But if you *really* want to,
186  * nothing is stopping you from creating a randutils::seed_seq_fe<624>.
187  *
188  * * As a consequence of the above properties, if all parts of the provided
189  * seed data are kept constant except one, and the remaining part is varied
190  * through K different states, K different output sequences will be
191  * produced.
192  *
193  * * Also, because the amount of entropy stored is fixed, this class never
194  * performs dynamic allocation and is free of the possibility of generating
195  * an exception.
196  *
197  * Ideas used to implement this code include hashing, a simple PCG generator
198  * based on an MCG base with an XorShift output function and permutation
199  * functions on tuples.
200  *
201  * More detail at
202  * http://www.pcg-random.org/posts/developing-a-seed_seq-alternative.html
203  */
204 
205  template<std::size_t count, typename IntRep = std::uint32_t>
206  struct seed_seq_fe
207  {
208  public:
209  CPP_assert(unsigned_integral<IntRep>);
210  typedef IntRep result_type;
211 
212  private:
213  static constexpr std::size_t mix_rounds = 1 + (count <= 2);
214 
215  static constexpr std::uint32_t INIT_A = 0x43b0d7e5;
216  static constexpr std::uint32_t MULT_A = 0x931e8875;
217 
218  static constexpr std::uint32_t INIT_B = 0x8b51f9dd;
219  static constexpr std::uint32_t MULT_B = 0x58f38ded;
220 
221  static constexpr std::uint32_t MIX_MULT_L = 0xca01f9dd;
222  static constexpr std::uint32_t MIX_MULT_R = 0x4973f715;
223  static constexpr std::uint32_t XSHIFT = sizeof(IntRep) * 8 / 2;
224 
225  std::array<IntRep, count> mixer_;
226 
227  template(typename I, typename S)(
228  requires input_iterator<I> AND sentinel_for<S, I> AND
229  convertible_to<iter_reference_t<I>, IntRep>)
230  void mix_entropy(I first, S last)
231  {
232  auto hash_const = INIT_A;
233  auto hash = [&](IntRep value) RANGES_INTENDED_MODULAR_ARITHMETIC {
234  value ^= hash_const;
235  hash_const *= MULT_A;
236  value *= hash_const;
237  value ^= value >> XSHIFT;
238  return value;
239  };
240  auto mix = [](IntRep x, IntRep y) RANGES_INTENDED_MODULAR_ARITHMETIC {
241  IntRep result = MIX_MULT_L * x - MIX_MULT_R * y;
242  result ^= result >> XSHIFT;
243  return result;
244  };
245 
246  for(auto & elem : mixer_)
247  {
248  if(first != last)
249  {
250  elem = hash(static_cast<IntRep>(*first));
251  ++first;
252  }
253  else
254  elem = hash(IntRep{0});
255  }
256  for(auto & src : mixer_)
257  for(auto & dest : mixer_)
258  if(&src != &dest)
259  dest = mix(dest, hash(src));
260  for(; first != last; ++first)
261  for(auto & dest : mixer_)
262  dest = mix(dest, hash(static_cast<IntRep>(*first)));
263  }
264 
265  public:
266  seed_seq_fe(const seed_seq_fe &) = delete;
267  void operator=(const seed_seq_fe &) = delete;
268 
269  template(typename T)(
270  requires convertible_to<T const &, IntRep>)
271  seed_seq_fe(std::initializer_list<T> init)
272  {
273  seed(init.begin(), init.end());
274  }
275 
276  template(typename I, typename S)(
277  requires input_iterator<I> AND sentinel_for<S, I> AND
278  convertible_to<iter_reference_t<I>, IntRep>)
279  seed_seq_fe(I first, S last)
280  {
281  seed(first, last);
282  }
283 
284  // generating functions
285  template(typename I, typename S)(
286  requires random_access_iterator<I> AND sentinel_for<S, I>)
287  RANGES_INTENDED_MODULAR_ARITHMETIC //
288  void generate(I first, S const last) const
289  {
290  auto src_begin = mixer_.begin();
291  auto src_end = mixer_.end();
292  auto src = src_begin;
293  auto hash_const = INIT_B;
294  for(; first != last; ++first)
295  {
296  auto dataval = *src;
297  if(++src == src_end)
298  src = src_begin;
299  dataval ^= hash_const;
300  hash_const *= MULT_B;
301  dataval *= hash_const;
302  dataval ^= dataval >> XSHIFT;
303  *first = dataval;
304  }
305  }
306 
307  constexpr std::size_t size() const
308  {
309  return count;
310  }
311 
312  template(typename O)(
313  requires weakly_incrementable<O> AND
314  indirectly_copyable<decltype(mixer_.begin()), O>)
315  RANGES_INTENDED_MODULAR_ARITHMETIC void param(O dest) const
316  {
317  constexpr IntRep INV_A = randutils::fast_exp(MULT_A, IntRep(-1));
318  constexpr IntRep MIX_INV_L =
319  randutils::fast_exp(MIX_MULT_L, IntRep(-1));
320 
321  auto mixer_copy = mixer_;
322  for(std::size_t round = 0; round < mix_rounds; ++round)
323  {
324  // Advance to the final value. We'll backtrack from that.
325  auto hash_const =
326  INIT_A * randutils::fast_exp(MULT_A, IntRep(count * count));
327 
328  for(auto src = mixer_copy.rbegin(); src != mixer_copy.rend();
329  ++src)
330  for(auto rdest = mixer_copy.rbegin();
331  rdest != mixer_copy.rend();
332  ++rdest)
333  if(src != rdest)
334  {
335  IntRep revhashed = *src;
336  auto mult_const = hash_const;
337  hash_const *= INV_A;
338  revhashed ^= hash_const;
339  revhashed *= mult_const;
340  revhashed ^= revhashed >> XSHIFT;
341  IntRep unmixed = *rdest;
342  unmixed ^= unmixed >> XSHIFT;
343  unmixed += MIX_MULT_R * revhashed;
344  unmixed *= MIX_INV_L;
345  *rdest = unmixed;
346  }
347 
348  for(auto i = mixer_copy.rbegin(); i != mixer_copy.rend(); ++i)
349  {
350  IntRep unhashed = *i;
351  unhashed ^= unhashed >> XSHIFT;
352  unhashed *= randutils::fast_exp(hash_const, IntRep(-1));
353  hash_const *= INV_A;
354  unhashed ^= hash_const;
355  *i = unhashed;
356  }
357  }
358  ranges::copy(mixer_copy, dest);
359  }
360 
361  template(typename I, typename S)(
362  requires input_iterator<I> AND sentinel_for<S, I> AND
363  convertible_to<iter_reference_t<I>, IntRep>)
364  void seed(I first, S last)
365  {
366  mix_entropy(first, last);
367  // For very small sizes, we do some additional mixing. For normal
368  // sizes, this loop never performs any iterations.
369  for(std::size_t i = 1; i < mix_rounds; ++i)
370  stir();
371  }
372 
373  seed_seq_fe & stir()
374  {
375  mix_entropy(mixer_.begin(), mixer_.end());
376  return *this;
377  }
378  };
379 
380  using seed_seq_fe128 = seed_seq_fe<4, std::uint32_t>;
381  using seed_seq_fe256 = seed_seq_fe<8, std::uint32_t>;
382 
384  //
385  // auto_seeded
386  //
388 
389  /*
390  * randutils::auto_seeded
391  *
392  * Extends a seed sequence class with a nondeterministic default
393  * constructor. Uses a variety of local sources of entropy to portably
394  * initialize any seed sequence to a good default state.
395  *
396  * In normal use, it's accessed via one of the following type aliases, which
397  * use seed_seq_fe128 and seed_seq_fe256 above.
398  *
399  * randutils::auto_seed_128
400  * randutils::auto_seed_256
401  *
402  * It's discussed in detail at
403  * http://www.pcg-random.org/posts/simple-portable-cpp-seed-entropy.html
404  * and its motivation (why you can't just use std::random_device) here
405  * http://www.pcg-random.org/posts/cpps-random_device.html
406  */
407 
408  template<typename SeedSeq>
409  struct auto_seeded : public SeedSeq
410  {
411  auto_seeded()
412  : auto_seeded(randutils::get_entropy())
413  {}
414  template<std::size_t N>
415  auto_seeded(std::array<std::uint32_t, N> const & seeds)
416  : SeedSeq(seeds.begin(), seeds.end())
417  {}
418  using SeedSeq::SeedSeq;
419 
420  const SeedSeq & base() const
421  {
422  return *this;
423  }
424  SeedSeq & base()
425  {
426  return *this;
427  }
428  };
429 
430  using auto_seed_128 = auto_seeded<seed_seq_fe128>;
431  using auto_seed_256 = auto_seeded<seed_seq_fe256>;
432  } // namespace randutils
433 
434  using default_URNG = meta::if_c<(sizeof(void *) >= sizeof(long long)),
435  std::mt19937_64, std::mt19937>;
436 
437 #if !RANGES_CXX_THREAD_LOCAL
438  template<typename URNG>
439  class sync_URNG : private URNG
440  {
441  mutable std::mutex mtx_;
442 
443  public:
444  using URNG::URNG;
445  sync_URNG() = default;
446  using typename URNG::result_type;
447  result_type operator()()
448  {
449  std::lock_guard<std::mutex> guard{mtx_};
450  return static_cast<URNG &>(*this)();
451  }
452  using URNG::max;
453  using URNG::min;
454  };
455  using default_random_engine = sync_URNG<default_URNG>;
456 #else
457  using default_random_engine = default_URNG;
458 #endif
459 
460  template<typename T = void>
461  default_random_engine & get_random_engine()
462  {
463  using Seeder = meta::if_c<(sizeof(default_URNG) > 16),
464  randutils::auto_seed_256,
465  randutils::auto_seed_128>;
466 
467 #if RANGES_CXX_THREAD_LOCAL >= RANGES_CXX_THREAD_LOCAL_11
468  static thread_local default_random_engine engine{Seeder{}.base()};
469 
470 #elif RANGES_CXX_THREAD_LOCAL
471  static __thread bool initialized = false;
472  static __thread meta::_t<std::aligned_storage<sizeof(default_random_engine),
473  alignof(default_random_engine)>>
474  storage;
475 
476  if(!initialized)
477  {
478  ::new(static_cast<void *>(&storage))
479  default_random_engine{Seeder{}.base()};
480  initialized = true;
481  }
482  auto & engine = reinterpret_cast<default_random_engine &>(storage);
483 #else
484  static default_random_engine engine{Seeder{}.base()};
485 #endif // RANGES_CXX_THREAD_LOCAL
486 
487  return engine;
488  }
489  } // namespace detail
491 } // namespace ranges
492 
493 RANGES_DIAGNOSTIC_POP
494 
495 #include <range/v3/detail/epilogue.hpp>
496 
497 #endif
CPP_concept indirectly_copyable
\concept indirectly_copyable
Definition: concepts.hpp:782
CPP_concept uniform_random_bit_generator
\concept uniform_random_bit_generator
Definition: random.hpp:99
CPP_requires(uniform_random_bit_generator_, requires()(Gen::min(), Gen::max()))
\concept uniform_random_bit_generator_
std::integral_constant< std::size_t, N > size_t
An integral constant wrapper for std::size_t.
Definition: meta.hpp:163
typename T::type _t
Type alias for T::type.
Definition: meta.hpp:141
front< Pair > first
Retrieve the first element of the pair Pair.
Definition: meta.hpp:2251
meta::size_t< L::size()> size
An integral constant wrapper that is the size of the meta::list L.
Definition: meta.hpp:1696
_t< detail::count_< L, T > > count
Count the number of times a type T appears in the list L.
Definition: meta.hpp:2725
Tiny meta-programming library.
std::size_t hash(const BasicJsonType &j)
hash a JSON value
Definition: hash.hpp:34