Line data Source code
1 : /****************************************************************************/
2 : // Eclipse SUMO, Simulation of Urban MObility; see https://eclipse.dev/sumo
3 : // Copyright (C) 2005-2025 German Aerospace Center (DLR) and others.
4 : // This program and the accompanying materials are made available under the
5 : // terms of the Eclipse Public License 2.0 which is available at
6 : // https://www.eclipse.org/legal/epl-2.0/
7 : // This Source Code may also be made available under the following Secondary
8 : // Licenses when the conditions for such availability set forth in the Eclipse
9 : // Public License 2.0 are satisfied: GNU General Public License, version 2
10 : // or later which is available at
11 : // https://www.gnu.org/licenses/old-licenses/gpl-2.0-standalone.html
12 : // SPDX-License-Identifier: EPL-2.0 OR GPL-2.0-or-later
13 : /****************************************************************************/
14 : /// @file RandHelper.h
15 : /// @author Daniel Krajzewicz
16 : /// @author Michael Behrisch
17 : /// @author Jakob Erdmann
18 : /// @date Fri, 29.04.2005
19 : ///
20 : //
21 : /****************************************************************************/
22 : #pragma once
23 : #include <config.h>
24 :
25 : #include <cassert>
26 : #include <cstring>
27 : #include <vector>
28 : #include <map>
29 : #include <random>
30 : #include <sstream>
31 : #include <iostream>
32 : #include <algorithm>
33 :
34 : // TODO make this configurable
35 : #define SAVE_ONLY_COUNT 1000000
36 :
37 : // ===========================================================================
38 : // helper function
39 : // ===========================================================================
40 :
41 : #ifdef __clang__
42 : __attribute__((no_sanitize("unsigned-integer-overflow")))
43 : #endif
44 : inline uint64_t splitmix64(const uint64_t seed) {
45 65106 : uint64_t z = (seed + 0x9e3779b97f4a7c15);
46 65106 : z = (z ^ (z >> 30)) * 0xbf58476d1ce4e5b9;
47 65106 : z = (z ^ (z >> 27)) * 0x94d049bb133111eb;
48 65106 : return z ^ (z >> 31);
49 : }
50 :
51 :
52 : // ===========================================================================
53 : // class declaration
54 : // ===========================================================================
55 :
56 : class OptionsCont;
57 :
58 : // ===========================================================================
59 : // class definitions
60 : // ===========================================================================
61 : /**
62 : * @class XoShiRo256PlusPlus
63 : * @brief A random number generator as proposed here https://prng.di.unimi.it/xoshiro256plusplus.c
64 : */
65 : class XoShiRo256PlusPlus {
66 : public:
67 : inline void seed(const uint64_t seed) {
68 : this->state[0] = splitmix64(splitmix64(seed));
69 : this->state[1] = splitmix64(this->state[0]);
70 : this->state[2] = splitmix64(this->state[1]);
71 : this->state[3] = splitmix64(this->state[2]);
72 : }
73 :
74 : inline uint64_t operator()() {
75 : const uint64_t result = rotl64(this->state[0] + this->state[3], 23) + this->state[0];
76 : const uint64_t t = this->state[1] << 17;
77 : this->state[2] ^= this->state[0];
78 : this->state[3] ^= this->state[1];
79 : this->state[1] ^= this->state[2];
80 : this->state[0] ^= this->state[3];
81 : this->state[2] ^= t;
82 : this->state[3] = rotl64(this->state[3], 45);
83 : return result;
84 : }
85 :
86 : void discard(unsigned long long skip) {
87 : while (skip-- > 0) {
88 : (*this)();
89 : }
90 : }
91 :
92 : friend std::ostream& operator<<(std::ostream& os, const XoShiRo256PlusPlus& r) {
93 : os << r.state[0] << r.state[1] << r.state[2] << r.state[3];
94 : return os;
95 : }
96 :
97 : friend std::istream& operator>>(std::istream& is, XoShiRo256PlusPlus& r) {
98 : is >> r.state[0] >> r.state[1] >> r.state[2] >> r.state[3];
99 : return is;
100 : }
101 :
102 :
103 : private:
104 : static inline uint64_t rotl64(const uint64_t x, const int k) {
105 : return (x << k) | (x >> (64 - k));
106 : }
107 :
108 : uint64_t state[4];
109 :
110 : };
111 :
112 :
113 : //class SumoRNG : public XoShiRo256PlusPlus {
114 7564051 : class SumoRNG : public std::mt19937 {
115 : public:
116 5702460 : SumoRNG(const std::string& _id) : id(_id) {}
117 :
118 : void setSeed(int _seed) {
119 2785718 : origSeed = _seed;
120 2785718 : seed(_seed);
121 2785718 : }
122 :
123 : unsigned long long int count = 0;
124 : int origSeed = default_seed;
125 : std::string id;
126 : };
127 :
128 :
129 : /**
130 : * @class RandHelper
131 : * @brief Utility functions for using a global, resetable random number generator
132 : */
133 : class RandHelper {
134 :
135 : public:
136 : /// @brief Initialises the given options container with random number options
137 : static void insertRandOptions(OptionsCont& oc);
138 :
139 : /// @brief Initialises the random number generator with hardware randomness or seed
140 : static void initRand(SumoRNG* which = nullptr, const bool random = false, const int seed = 23423);
141 :
142 : /// @brief Reads the given random number options and initialises the random number generator in accordance
143 : static void initRandGlobal(SumoRNG* which = nullptr);
144 :
145 : static int getSeed(SumoRNG* which = nullptr);
146 :
147 : /// @brief Returns a random real number in [0, 1)
148 : static double rand(SumoRNG* rng = nullptr);
149 :
150 : /// @brief Returns a random real number in [0, maxV)
151 : static inline double rand(double maxV, SumoRNG* rng = nullptr) {
152 51472903 : return maxV * rand(rng);
153 : }
154 :
155 : /// @brief Returns a random real number in [minV, maxV)
156 : static inline double rand(double minV, double maxV, SumoRNG* rng = nullptr) {
157 194962 : return minV + (maxV - minV) * rand(rng);
158 : }
159 :
160 : /// @brief Returns a random integer in [0, maxV-1]
161 81857781 : static inline int rand(int maxV, SumoRNG* rng = nullptr) {
162 83219259 : if (rng == nullptr) {
163 : rng = &myRandomNumberGenerator;
164 : }
165 81857781 : unsigned int usedBits = maxV - 1;
166 81857781 : usedBits |= usedBits >> 1;
167 81857781 : usedBits |= usedBits >> 2;
168 81857781 : usedBits |= usedBits >> 4;
169 81857781 : usedBits |= usedBits >> 8;
170 81857781 : usedBits |= usedBits >> 16;
171 :
172 : // Draw numbers until one is found in [0, maxV-1]
173 : int result;
174 : do {
175 84220498 : result = (*rng)() & usedBits;
176 84220498 : rng->count++;
177 82798753 : } while (result >= maxV);
178 81857781 : return result;
179 : }
180 :
181 : /// @brief Returns a random integer in [minV, maxV-1]
182 : static inline int rand(int minV, int maxV, SumoRNG* rng = nullptr) {
183 92735 : return minV + rand(maxV - minV, rng);
184 : }
185 :
186 : /// @brief Returns a random 64 bit integer in [0, maxV-1]
187 14252 : static inline long long int rand(long long int maxV, SumoRNG* rng = nullptr) {
188 14252 : if (maxV <= std::numeric_limits<int>::max()) {
189 14252 : return rand((int)maxV, rng);
190 : }
191 0 : if (rng == nullptr) {
192 : rng = &myRandomNumberGenerator;
193 : }
194 0 : unsigned long long int usedBits = maxV - 1;
195 0 : usedBits |= usedBits >> 1;
196 0 : usedBits |= usedBits >> 2;
197 0 : usedBits |= usedBits >> 4;
198 0 : usedBits |= usedBits >> 8;
199 0 : usedBits |= usedBits >> 16;
200 0 : usedBits |= usedBits >> 32;
201 :
202 : // Draw numbers until one is found in [0, maxV-1]
203 : long long int result;
204 : do {
205 0 : result = (((unsigned long long int)(*rng)() << 32) | (*rng)()) & usedBits; // toss unused bits to shorten search
206 0 : rng->count += 2;
207 0 : } while (result >= maxV);
208 : return result;
209 : }
210 :
211 : /// @brief Returns a random 64 bit integer in [minV, maxV-1]
212 : static inline long long int rand(long long int minV, long long int maxV, SumoRNG* rng = nullptr) {
213 14032 : return minV + rand(maxV - minV, rng);
214 : }
215 :
216 : /// @brief Access to a random number from a normal distribution
217 : static double randNorm(double mean, double variance, SumoRNG* rng = nullptr);
218 :
219 : /// @brief Access to a random number from an exponential distribution
220 : static double randExp(double rate, SumoRNG* rng = nullptr);
221 :
222 : /// @brief Returns a random element from the given vector
223 : template<class T>
224 : static inline const T&
225 : getRandomFrom(const std::vector<T>& v, SumoRNG* rng = nullptr) {
226 : assert(v.size() > 0);
227 81670460 : return v[rand((int)v.size(), rng)];
228 : }
229 :
230 : /// @brief save rng state to string
231 781 : static std::string saveState(SumoRNG* rng = nullptr) {
232 781 : if (rng == nullptr) {
233 : rng = &myRandomNumberGenerator;
234 : }
235 781 : std::ostringstream oss;
236 781 : oss << rng->count;
237 781 : if (rng->count >= SAVE_ONLY_COUNT) {
238 1 : oss << " " << (*rng);
239 : }
240 781 : return oss.str();
241 781 : }
242 :
243 : /// @brief load rng state from string
244 781 : static void loadState(const std::string& state, SumoRNG* rng = nullptr) {
245 781 : if (rng == nullptr) {
246 : rng = &myRandomNumberGenerator;
247 : }
248 781 : std::istringstream iss(state);
249 781 : iss >> rng->count;
250 781 : if (rng->count < SAVE_ONLY_COUNT) {
251 780 : rng->discard(rng->count);
252 : } else {
253 1 : iss >> (*rng);
254 : }
255 781 : }
256 :
257 : template<class T>
258 0 : static void shuffle(std::vector<T>& v, SumoRNG* rng = nullptr) {
259 0 : for (int i = (int)(v.size() - 1); i > 0; --i) {
260 0 : std::swap(*(v.begin() + i), *(v.begin() + rand(i, rng)));
261 : }
262 0 : }
263 :
264 : static long long int count() {
265 : return myRandomNumberGenerator.count;
266 : }
267 :
268 : /// @brief return a value scrambled value from [0, 1]
269 : static double randHash(long long int x) {
270 65106 : const uint64_t h = splitmix64(x) >> 11;
271 65106 : return (double)h / (double(1ULL << 53));
272 : }
273 :
274 : /// @brief string hashing adapted from https://en.wikipedia.org/wiki/MurmurHash
275 : #ifdef __clang__
276 : __attribute__((no_sanitize("integer"))) // left-shift and unsigned-integer-overflow
277 : #endif
278 6196637 : static uint32_t murmur3_32(const std::string& key2, int seed) {
279 : const uint8_t* key = reinterpret_cast<const uint8_t*>(key2.data());
280 : const size_t len = key2.size();
281 6196637 : uint32_t h = seed;
282 : uint32_t k;
283 : /* Read in groups of 4. */
284 15773889 : for (size_t i = len >> 2; i; i--) {
285 : memcpy(&k, key, sizeof(uint32_t));
286 9577252 : key += sizeof(uint32_t);
287 9577252 : h ^= murmur_32_scramble(k);
288 9577252 : h = (h << 13) | (h >> 19);
289 9577252 : h = h * 5 + 0xe6546b64;
290 : }
291 : /* Read the rest. */
292 6196637 : k = 0;
293 16043136 : for (size_t i = len & 3; i; i--) {
294 9846499 : k <<= 8;
295 9846499 : k |= key[i - 1];
296 : }
297 6196637 : h ^= murmur_32_scramble(k);
298 : /* Finalize. */
299 6196637 : h ^= len;
300 6196637 : h ^= h >> 16;
301 6196637 : h *= 0x85ebca6b;
302 6196637 : h ^= h >> 13;
303 6196637 : h *= 0xc2b2ae35;
304 6196637 : h ^= h >> 16;
305 6196637 : return h;
306 : }
307 :
308 :
309 : protected:
310 : /// @brief the default random number generator to use
311 : static SumoRNG myRandomNumberGenerator;
312 :
313 : /// @brief helper function for murmur_32_scramble from https://en.wikipedia.org/wiki/MurmurHash
314 : #ifdef __clang__
315 : __attribute__((no_sanitize("integer"))) // left-shift and unsigned-integer-overflow
316 : #endif
317 : static inline uint32_t murmur_32_scramble(uint32_t k) {
318 15773889 : k *= 0xcc9e2d51;
319 15773889 : k = (k << 15) | (k >> 17);
320 9577252 : k *= 0x1b873593;
321 : return k;
322 : }
323 :
324 : };
|
