rng.cpp
Go to the documentation of this file.
1 
9 #include "rng.h"
10 #include <argos3/core/utility/configuration/argos_exception.h>
11 #include <argos3/core/utility/logging/argos_log.h>
12 #include <cstring>
13 #include <limits>
14 #include <cmath>
15 
16 namespace argos {
17 
18  /****************************************/
19  /****************************************/
20 
21  /* Period parameters */
22  static const SInt32 N = 624;
23  static const SInt32 M = 397;
24  static const UInt32 MATRIX_A = 0x9908b0dfUL; /* constant vector a */
25  static const UInt32 UPPER_MASK = 0x80000000UL; /* most significant w-r bits */
26  static const UInt32 LOWER_MASK = 0x7fffffffUL; /* least significant r bits */
27  static const CRange<UInt32> INT_RANGE = CRange<UInt32>(0, 0xFFFFFFFFUL);
28 
29  std::map<std::string, CRandom::CCategory*> CRandom::m_mapCategories;
30 
31  /* Checks that a category exists. It internally creates an iterator that points to the category, if found. */
32 #define CHECK_CATEGORY(category) \
33  std::map<std::string, CCategory*>::iterator itCategory = m_mapCategories.find(category); \
34  if(itCategory == m_mapCategories.end()) { \
35  THROW_ARGOSEXCEPTION("CRandom:: can't find category \"" << category << "\"."); \
36  }
37 
38  /****************************************/
39  /****************************************/
40 
41  CRandom::CRNG::CRNG(UInt32 un_seed) :
42  m_unSeed(un_seed),
43  m_punState(new UInt32[N]),
44  m_nIndex(N+1) {
45  Reset();
46  }
47 
48  /****************************************/
49  /****************************************/
50 
51  CRandom::CRNG::CRNG(const CRNG& c_rng) :
52  m_unSeed(c_rng.m_unSeed),
53  m_punState(new UInt32[N]),
54  m_nIndex(c_rng.m_nIndex) {
55  ::memcpy(m_punState, c_rng.m_punState, N * sizeof(UInt32));
56  }
57 
58  /****************************************/
59  /****************************************/
60 
61  CRandom::CRNG::~CRNG() {
62  delete[] m_punState;
63  }
64 
65  /****************************************/
66  /****************************************/
67 
68  void CRandom::CRNG::Reset() {
69  m_punState[0]= m_unSeed & 0xffffffffUL;
70  for (m_nIndex = 1; m_nIndex < N; ++m_nIndex) {
71  m_punState[m_nIndex] =
72  (1812433253UL * (m_punState[m_nIndex-1] ^ (m_punState[m_nIndex-1] >> 30)) + m_nIndex);
73  m_punState[m_nIndex] &= 0xffffffffUL;
74  }
75  }
76 
77  /****************************************/
78  /****************************************/
79 
80  bool CRandom::CRNG::Bernoulli(Real f_true) {
81  return Uniform32bit() < f_true * INT_RANGE.GetMax();
82  }
83 
84  /****************************************/
85  /****************************************/
86 
87  CRadians CRandom::CRNG::Uniform(const CRange<CRadians>& c_range) {
88  CRadians cRetVal;
89  INT_RANGE.MapValueIntoRange(cRetVal, Uniform32bit(), c_range);
90  return cRetVal;
91  }
92 
93  /****************************************/
94  /****************************************/
95 
96  Real CRandom::CRNG::Uniform(const CRange<Real>& c_range) {
97  Real fRetVal;
98  INT_RANGE.MapValueIntoRange(fRetVal, Uniform32bit(), c_range);
99  return fRetVal;
100  }
101 
102  /****************************************/
103  /****************************************/
104 
105  SInt32 CRandom::CRNG::Uniform(const CRange<SInt32>& c_range) {
106  SInt32 nRetVal;
107  INT_RANGE.MapValueIntoRange(nRetVal, Uniform32bit(), c_range);
108  return nRetVal;
109  }
110 
111  /****************************************/
112  /****************************************/
113 
114  UInt32 CRandom::CRNG::Uniform(const CRange<UInt32>& c_range) {
115  UInt32 unRetVal;
116  INT_RANGE.MapValueIntoRange(unRetVal, Uniform32bit(), c_range);
117  return unRetVal;
118  }
119 
120  /****************************************/
121  /****************************************/
122 
123  Real CRandom::CRNG::Exponential(Real f_mean) {
124  static CRange<Real> fRange(0.0f, 1.0f);
125  return -Log(Uniform(fRange)) * f_mean;
126  }
127 
128  /****************************************/
129  /****************************************/
130 
131  UInt32 CRandom::CRNG::Poisson(Real f_mean)
132  {
133  /* This is the Devroye method.
134  see http://www.nrbook.com/devroye/Devroye_files/chapter_six.pdf
135  */
136  UInt32 unRetVal = 0;
137  Real fP, fS, fU;
138  fP = Exp(-f_mean);
139  fS = fP;
140  fU = Uniform(CRange<Real>(0.0f, 1.0f));
141  while(fU > fS) {
142  unRetVal += 1;
143  fP *= f_mean / unRetVal;
144  fS += fP;
145  }
146  return unRetVal;
147  }
148 
149  /****************************************/
150  /****************************************/
151 
152  Real CRandom::CRNG::Gaussian(Real f_std_dev,
153  Real f_mean) {
154  /* This is the Box-Muller method in its cartesian variant
155  see http://www.dspguru.com/dsp/howtos/how-to-generate-white-gaussian-noise
156  */
157  static CRange<Real> fRange(-1.0f, 1.0f);
158  Real fNum1, fNum2;
159  Real fSquare;
160  do {
161  fNum1 = Uniform(fRange);
162  fNum2 = Uniform(fRange);
163  fSquare = fNum1 * fNum1 + fNum2 * fNum2;
164  } while(fSquare >= 1);
165  return f_mean + f_std_dev * fNum1 * Sqrt(-2.0f * Log(fSquare) / fSquare);
166  }
167 
168  /****************************************/
169  /****************************************/
170 
171  Real CRandom::CRNG::Rayleigh(Real f_sigma) {
172  /* Draw a number uniformly from (0,1) --- bounds excluded */
173  static CRange<Real> cUnitRange(0.0f, 1.0f);
174  Real fValue;
175  do {
176  fValue = Uniform(cUnitRange);
177  }
178  while(! cUnitRange.WithinMinBoundExcludedMaxBoundExcluded(fValue));
179  /* Calculate the value to return from the definition of Rayleigh distribution
180  * http://en.wikipedia.org/wiki/Rayleigh_distribution#Generating_Rayleigh-distributed_random_variates
181  */
182  return f_sigma * Sqrt(-2.0f * Log(fValue));
183  }
184 
185 /****************************************/
186  /****************************************/
187 
188  Real CRandom::CRNG::Lognormal(Real f_sigma, Real f_mu) {
189  /* Draw a number uniformly from (0,1) */
190  Real fValue;
191  fValue = Gaussian(1,0);
192  /* Calculate the value to return from the definition of Lognormal distribution
193  * http://en.wikipedia.org/wiki/Log-normal_distribution#Generating_log-normally_distributed_random_variates
194  */
195  return std::exp(f_mu + f_sigma * fValue);
196  }
197 
198  /****************************************/
199  /****************************************/
200 
201  UInt32 CRandom::CRNG::Uniform32bit() {
202  UInt32 y;
203  static UInt32 mag01[2] = { 0x0UL, MATRIX_A };
204  /* mag01[x] = x * MATRIX_A for x=0,1 */
205 
206  if (m_nIndex >= N) { /* generate N words at one time */
207  SInt32 kk;
208  for (kk = 0; kk < N - M; ++kk) {
209  y = (m_punState[kk] & UPPER_MASK) | (m_punState[kk+1] & LOWER_MASK);
210  m_punState[kk] = m_punState[kk+M] ^ (y >> 1) ^ mag01[y & 0x1UL];
211  }
212  for (; kk < N - 1; ++kk) {
213  y = (m_punState[kk] & UPPER_MASK) | (m_punState[kk+1] & LOWER_MASK);
214  m_punState[kk] = m_punState[kk+(M-N)] ^ (y >> 1) ^ mag01[y & 0x1UL];
215  }
216  y = (m_punState[N-1] & UPPER_MASK) | (m_punState[0] & LOWER_MASK);
217  m_punState[N-1] = m_punState[M-1] ^ (y >> 1) ^ mag01[y & 0x1UL];
218 
219  m_nIndex = 0;
220  }
221 
222  y = m_punState[m_nIndex++];
223 
224  /* Tempering */
225  y ^= (y >> 11);
226  y ^= (y << 7) & 0x9d2c5680UL;
227  y ^= (y << 15) & 0xefc60000UL;
228  y ^= (y >> 18);
229 
230  return y;
231  }
232 
233  /****************************************/
234  /****************************************/
235 
236  CRandom::CCategory::CCategory(const std::string& str_id,
237  UInt32 un_seed) :
238  m_strId(str_id),
239  m_unSeed(un_seed),
240  m_cSeeder(un_seed),
241  m_cSeedRange(1, std::numeric_limits<UInt32>::max()) {}
242 
243  /****************************************/
244  /****************************************/
245 
246  CRandom::CCategory::~CCategory() {
247  while(! m_vecRNGList.empty()) {
248  delete m_vecRNGList.back();
249  m_vecRNGList.pop_back();
250  }
251  }
252 
253  /****************************************/
254  /****************************************/
255 
256  void CRandom::CCategory::SetSeed(UInt32 un_seed) {
257  m_unSeed = un_seed;
258  m_cSeeder.SetSeed(m_unSeed);
259  }
260 
261  /****************************************/
262  /****************************************/
263 
264  CRandom::CRNG* CRandom::CCategory::CreateRNG() {
265  /* Get seed from internal RNG */
266  UInt32 unSeed = m_cSeeder.Uniform(m_cSeedRange);
267  /* Create new RNG */
268  m_vecRNGList.push_back(new CRNG(unSeed));
269  return m_vecRNGList.back();
270  }
271 
272  /****************************************/
273  /****************************************/
274 
275  void CRandom::CCategory::ResetRNGs() {
276  /* Reset internal RNG */
277  m_cSeeder.Reset();
278  ReseedRNGs();
279  /* Reset the RNGs */
280  for(size_t i = 0; i < m_vecRNGList.size(); ++i) {
281  m_vecRNGList[i]->Reset();
282  }
283  }
284 
285  /****************************************/
286  /****************************************/
287 
288  void CRandom::CCategory::ReseedRNGs() {
289  for(size_t i = 0; i < m_vecRNGList.size(); ++i) {
290  /* Get seed from internal RNG */
291  m_vecRNGList[i]->SetSeed(m_cSeeder.Uniform(m_cSeedRange));
292  }
293  }
294 
295  /****************************************/
296  /****************************************/
297 
298  bool CRandom::CreateCategory(const std::string& str_category,
299  UInt32 un_seed) {
300  /* Is there a category already? */
301  auto itCategory = m_mapCategories.find(str_category);
302  if(itCategory == m_mapCategories.end()) {
303  /* No, create it */
304  m_mapCategories.insert(
305  std::pair<std::string,
306  CRandom::CCategory*>(str_category,
307  new CRandom::CCategory(str_category,
308  un_seed)));
309  return true;
310  }
311  return false;
312  }
313 
314  /****************************************/
315  /****************************************/
316 
317  CRandom::CCategory& CRandom::GetCategory(const std::string& str_category) {
318  CHECK_CATEGORY(str_category);
319  return *(itCategory->second);
320  }
321 
322  /****************************************/
323  /****************************************/
324 
325  bool CRandom::ExistsCategory(const std::string& str_category) {
326  try {
327  CHECK_CATEGORY(str_category);
328  return true;
329  }
330  catch(CARGoSException& ex) {
331  return false;
332  }
333  }
334 
335  /****************************************/
336  /****************************************/
337 
338  void CRandom::RemoveCategory(const std::string& str_category) {
339  CHECK_CATEGORY(str_category);
340  delete itCategory->second;
341  m_mapCategories.erase(itCategory);
342  }
343 
344  /****************************************/
345  /****************************************/
346 
347  CRandom::CRNG* CRandom::CreateRNG(const std::string& str_category) {
348  CHECK_CATEGORY(str_category);
349  return itCategory->second->CreateRNG();
350  }
351 
352  /****************************************/
353  /****************************************/
354 
355  UInt32 CRandom::GetSeedOf(const std::string& str_category) {
356  CHECK_CATEGORY(str_category);
357  return itCategory->second->GetSeed();
358  }
359 
360  /****************************************/
361  /****************************************/
362 
363  void CRandom::SetSeedOf(const std::string& str_category,
364  UInt32 un_seed) {
365  CHECK_CATEGORY(str_category);
366  itCategory->second->SetSeed(un_seed);
367  }
368 
369  /****************************************/
370  /****************************************/
371 
372  void CRandom::Reset() {
373  for(auto itCategory = m_mapCategories.begin();
374  itCategory != m_mapCategories.end();
375  ++itCategory) {
376  itCategory->second->ResetRNGs();
377  }
378  }
379 
380  /****************************************/
381  /****************************************/
382 
383 }
SInt32
signed int SInt32
32-bit signed integer.
Definition: datatypes.h:93
UInt32
unsigned int UInt32
32-bit unsigned integer.
Definition: datatypes.h:97
Real
float Real
Collects all ARGoS code.
Definition: datatypes.h:39
Exp
#define Exp
Definition: general.h:65
Sqrt
#define Sqrt
Definition: general.h:64
Log
#define Log
Definition: general.h:63
CHECK_CATEGORY
#define CHECK_CATEGORY(category)
Definition: rng.cpp:32
M
#define M(row, col)
argos
The namespace containing all the ARGoS related code.
Definition: ci_actuator.h:12
argos::CARGoSException
The exception that wraps all errors in ARGoS.
Definition: argos_exception.h:61
argos::CRadians
It defines the basic type CRadians, used to store an angle value in radians.
Definition: angles.h:42
argos::CRange::WithinMinBoundExcludedMaxBoundExcluded
bool WithinMinBoundExcludedMaxBoundExcluded(const T &t_value) const
Definition: range.h:93
argos::CRange::GetMax
T GetMax() const
Definition: range.h:48
argos::CRange::MapValueIntoRange
void MapValueIntoRange(U &t_output_value, const T &t_input_value, const CRange< U > &c_range) const
Definition: range.h:109
argos::CRandom::GetSeedOf
static UInt32 GetSeedOf(const std::string &str_category)
Returns the seed of the wanted category.
Definition: rng.cpp:355
argos::CRandom::Reset
static void Reset()
Resets all the RNG categories.
Definition: rng.cpp:372
argos::CRandom::CreateRNG
static CRNG * CreateRNG(const std::string &str_category)
Creates a new RNG inside the given category.
Definition: rng.cpp:347
argos::CRandom::RemoveCategory
static void RemoveCategory(const std::string &str_category)
Removes the wanted category.
Definition: rng.cpp:338
argos::CRandom::SetSeedOf
static void SetSeedOf(const std::string &str_category, UInt32 un_seed)
Sets the new seed of the wanted category.
Definition: rng.cpp:363
argos::CRandom::CreateCategory
static bool CreateCategory(const std::string &str_category, UInt32 un_seed)
Creates a new category.
Definition: rng.cpp:298
argos::CRandom::GetCategory
static CCategory & GetCategory(const std::string &str_category)
Returns a reference to the wanted category.
Definition: rng.cpp:317
argos::CRandom::ExistsCategory
static bool ExistsCategory(const std::string &str_category)
Returns true if the given category exists in the pool.
Definition: rng.cpp:325
argos::CRandom::CRNG
The RNG.
Definition: rng.h:90
argos::CRandom::CRNG::Poisson
UInt32 Poisson(Real f_mean)
Returns a random value from a Poisson distribution.
Definition: rng.cpp:131
argos::CRandom::CRNG::~CRNG
virtual ~CRNG()
Class destructor.
Definition: rng.cpp:61
argos::CRandom::CRNG::Bernoulli
bool Bernoulli(Real f_true=0.5)
Returns a random value from a Bernoulli distribution.
Definition: rng.cpp:80
argos::CRandom::CRNG::CRNG
CRNG(UInt32 un_seed)
Class constructor.
Definition: rng.cpp:41
argos::CRandom::CRNG::Reset
void Reset()
Reset the RNG.
Definition: rng.cpp:68
argos::CRandom::CRNG::Lognormal
Real Lognormal(Real f_sigma, Real f_mu)
Returns a random value from a Lognormal distribution.
Definition: rng.cpp:188
argos::CRandom::CRNG::Rayleigh
Real Rayleigh(Real f_sigma)
Returns a random value from a Rayleigh distribution.
Definition: rng.cpp:171
argos::CRandom::CRNG::Exponential
Real Exponential(Real f_mean)
Returns a random value from an exponential distribution.
Definition: rng.cpp:123
argos::CRandom::CRNG::Gaussian
Real Gaussian(Real f_std_dev, Real f_mean=0.0f)
Returns a random value from a Gaussian distribution.
Definition: rng.cpp:152
argos::CRandom::CRNG::Uniform
CRadians Uniform(const CRange< CRadians > &c_range)
Returns a random value from a uniform distribution.
Definition: rng.cpp:87
argos::CRandom::CCategory
The RNG category.
Definition: rng.h:244
argos::CRandom::CCategory::ResetRNGs
void ResetRNGs()
Resets the RNGs in this category.
Definition: rng.cpp:275
argos::CRandom::CCategory::CCategory
CCategory(const std::string &str_id, UInt32 un_seed)
Class constructor.
Definition: rng.cpp:236
argos::CRandom::CCategory::CreateRNG
CRNG * CreateRNG()
Creates a new RNG inside this category.
Definition: rng.cpp:264
argos::CRandom::CCategory::~CCategory
virtual ~CCategory()
Class destructor.
Definition: rng.cpp:246
argos::CRandom::CCategory::SetSeed
void SetSeed(UInt32 un_seed)
Sets the new seed of the category.
Definition: rng.cpp:256
argos::CRandom::CCategory::ReseedRNGs
void ReseedRNGs()
Sets new seed for the RNGs in this category.
Definition: rng.cpp:288