Line data Source code
1 : /****************************************************************************/
2 : // Eclipse SUMO, Simulation of Urban MObility; see https://eclipse.dev/sumo
3 : // Copyright (C) 2005-2024 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 <vector>
27 : #include <map>
28 : #include <random>
29 : #include <sstream>
30 : #include <iostream>
31 : #include <algorithm>
32 :
33 :
34 : // ===========================================================================
35 : // class declaration
36 : // ===========================================================================
37 :
38 : class OptionsCont;
39 :
40 : // ===========================================================================
41 : // class definitions
42 : // ===========================================================================
43 : /**
44 : * @class XoShiRo256PlusPlus
45 : * @brief A random number generator as proposed here https://prng.di.unimi.it/xoshiro256plusplus.c
46 : */
47 : class XoShiRo256PlusPlus {
48 : public:
49 : inline void seed(const uint64_t seed) {
50 : this->state[0] = splitmix64(splitmix64(seed));
51 : this->state[1] = splitmix64(this->state[0]);
52 : this->state[2] = splitmix64(this->state[1]);
53 : this->state[3] = splitmix64(this->state[2]);
54 : }
55 :
56 : inline uint64_t operator()() {
57 : const uint64_t result = rotl64(this->state[0] + this->state[3], 23) + this->state[0];
58 : const uint64_t t = this->state[1] << 17;
59 : this->state[2] ^= this->state[0];
60 : this->state[3] ^= this->state[1];
61 : this->state[1] ^= this->state[2];
62 : this->state[0] ^= this->state[3];
63 : this->state[2] ^= t;
64 : this->state[3] = rotl64(this->state[3], 45);
65 : return result;
66 : }
67 :
68 : void discard(unsigned long long skip) {
69 : while (skip-- > 0) {
70 : (*this)();
71 : }
72 : }
73 :
74 : friend std::ostream& operator<<(std::ostream& os, const XoShiRo256PlusPlus& r) {
75 : os << r.state[0] << r.state[1] << r.state[2] << r.state[3];
76 : return os;
77 : }
78 :
79 : friend std::istream& operator>>(std::istream& is, XoShiRo256PlusPlus& r) {
80 : is >> r.state[0] >> r.state[1] >> r.state[2] >> r.state[3];
81 : return is;
82 : }
83 :
84 :
85 : private:
86 : static inline uint64_t rotl64(const uint64_t x, const int k) {
87 : return (x << k) | (x >> (64 - k));
88 : }
89 :
90 : static inline uint64_t splitmix64(const uint64_t seed) {
91 : uint64_t z = (seed + 0x9e3779b97f4a7c15);
92 : z = (z ^ (z >> 30)) * 0xbf58476d1ce4e5b9;
93 : z = (z ^ (z >> 27)) * 0x94d049bb133111eb;
94 : return z ^ (z >> 31);
95 : }
96 :
97 : uint64_t state[4];
98 :
99 : };
100 :
101 :
102 : //class SumoRNG : public XoShiRo256PlusPlus {
103 8295448 : class SumoRNG : public std::mt19937 {
104 : public:
105 6258094 : SumoRNG(const std::string& _id) : id(_id) {}
106 :
107 : unsigned long long int count = 0;
108 : std::string id;
109 : };
110 :
111 :
112 : /**
113 : * @class RandHelper
114 : * @brief Utility functions for using a global, resetable random number generator
115 : */
116 : class RandHelper {
117 :
118 : public:
119 : /// @brief Initialises the given options container with random number options
120 : static void insertRandOptions(OptionsCont& oc);
121 :
122 : /// @brief Initialises the random number generator with hardware randomness or seed
123 : static void initRand(SumoRNG* which = nullptr, const bool random = false, const int seed = 23423);
124 :
125 : /// @brief Reads the given random number options and initialises the random number generator in accordance
126 : static void initRandGlobal(SumoRNG* which = nullptr);
127 :
128 : /// @brief Returns a random real number in [0, 1)
129 : static double rand(SumoRNG* rng = nullptr);
130 :
131 : /// @brief Returns a random real number in [0, maxV)
132 : static inline double rand(double maxV, SumoRNG* rng = nullptr) {
133 50869085 : return maxV * rand(rng);
134 : }
135 :
136 : /// @brief Returns a random real number in [minV, maxV)
137 : static inline double rand(double minV, double maxV, SumoRNG* rng = nullptr) {
138 200376 : return minV + (maxV - minV) * rand(rng);
139 : }
140 :
141 : /// @brief Returns a random integer in [0, maxV-1]
142 80288197 : static inline int rand(int maxV, SumoRNG* rng = nullptr) {
143 81783850 : if (rng == nullptr) {
144 : rng = &myRandomNumberGenerator;
145 : }
146 80288197 : unsigned int usedBits = maxV - 1;
147 80288197 : usedBits |= usedBits >> 1;
148 80288197 : usedBits |= usedBits >> 2;
149 80288197 : usedBits |= usedBits >> 4;
150 80288197 : usedBits |= usedBits >> 8;
151 80288197 : usedBits |= usedBits >> 16;
152 :
153 : // Draw numbers until one is found in [0, maxV-1]
154 : int result;
155 : do {
156 82459522 : result = (*rng)() & usedBits;
157 82459522 : rng->count++;
158 80903602 : } while (result >= maxV);
159 80288197 : return result;
160 : }
161 :
162 : /// @brief Returns a random integer in [minV, maxV-1]
163 : static inline int rand(int minV, int maxV, SumoRNG* rng = nullptr) {
164 89099 : return minV + rand(maxV - minV, rng);
165 : }
166 :
167 : /// @brief Returns a random 64 bit integer in [0, maxV-1]
168 14252 : static inline long long int rand(long long int maxV, SumoRNG* rng = nullptr) {
169 14252 : if (maxV <= std::numeric_limits<int>::max()) {
170 14252 : return rand((int)maxV, rng);
171 : }
172 0 : if (rng == nullptr) {
173 : rng = &myRandomNumberGenerator;
174 : }
175 0 : unsigned long long int usedBits = maxV - 1;
176 0 : usedBits |= usedBits >> 1;
177 0 : usedBits |= usedBits >> 2;
178 0 : usedBits |= usedBits >> 4;
179 0 : usedBits |= usedBits >> 8;
180 0 : usedBits |= usedBits >> 16;
181 0 : usedBits |= usedBits >> 32;
182 :
183 : // Draw numbers until one is found in [0, maxV-1]
184 : long long int result;
185 : do {
186 0 : result = (((unsigned long long int)(*rng)() << 32) | (*rng)()) & usedBits; // toss unused bits to shorten search
187 0 : rng->count += 2;
188 0 : } while (result >= maxV);
189 : return result;
190 : }
191 :
192 : /// @brief Returns a random 64 bit integer in [minV, maxV-1]
193 : static inline long long int rand(long long int minV, long long int maxV, SumoRNG* rng = nullptr) {
194 14032 : return minV + rand(maxV - minV, rng);
195 : }
196 :
197 : /// @brief Access to a random number from a normal distribution
198 : static double randNorm(double mean, double variance, SumoRNG* rng = nullptr);
199 :
200 : /// @brief Access to a random number from an exponential distribution
201 : static double randExp(double rate, SumoRNG* rng = nullptr);
202 :
203 : /// @brief Returns a random element from the given vector
204 : template<class T>
205 : static inline const T&
206 : getRandomFrom(const std::vector<T>& v, SumoRNG* rng = nullptr) {
207 : assert(v.size() > 0);
208 80114836 : return v[rand((int)v.size(), rng)];
209 : }
210 :
211 : /// @brief save rng state to string
212 426 : static std::string saveState(SumoRNG* rng = nullptr) {
213 426 : if (rng == nullptr) {
214 : rng = &myRandomNumberGenerator;
215 : }
216 426 : std::ostringstream oss;
217 426 : if (rng->count < 1000000) { // TODO make this configurable
218 : oss << rng->count;
219 : } else {
220 0 : oss << (*rng);
221 : }
222 426 : return oss.str();
223 426 : }
224 :
225 : /// @brief load rng state from string
226 426 : static void loadState(const std::string& state, SumoRNG* rng = nullptr) {
227 426 : if (rng == nullptr) {
228 : rng = &myRandomNumberGenerator;
229 : }
230 426 : std::istringstream iss(state);
231 426 : if (state.size() < 10) {
232 426 : iss >> rng->count;
233 426 : rng->discard(rng->count);
234 : } else {
235 0 : iss >> (*rng);
236 : }
237 426 : }
238 :
239 : template<class T>
240 0 : static void shuffle(std::vector<T>& v, SumoRNG* rng = nullptr) {
241 0 : for (int i = (int)(v.size() - 1); i > 0; --i) {
242 0 : std::swap(*(v.begin() + i), *(v.begin() + rand(i, rng)));
243 : }
244 0 : }
245 :
246 : protected:
247 : /// @brief the default random number generator to use
248 : static SumoRNG myRandomNumberGenerator;
249 :
250 : };
|
