Arcane  v3.16.8.0
Documentation développeur
Chargement...
Recherche...
Aucune correspondance
MersenneTwister.h
1// -*- tab-width: 2; indent-tabs-mode: nil; coding: utf-8-with-signature -*-
2//-----------------------------------------------------------------------------
3// Copyright 2000-2025 CEA (www.cea.fr) IFPEN (www.ifpenergiesnouvelles.com)
4// See the top-level COPYRIGHT file for details.
5// SPDX-License-Identifier: Apache-2.0
6//-----------------------------------------------------------------------------
7/*---------------------------------------------------------------------------*/
8/* MersenneTwister.h (C) 2000-2025 */
9/* */
10/* Ce fichier definit le patron de classe MersenneTwister ainsi que deux */
11/* classes associees mt19937 et mt11213b. Il est une version adapte a TROLL */
12/* du fichier MersenneTwister.hpp provenant de la bibliotheque BOOST */
13/*---------------------------------------------------------------------------*/
14/*---------------------------------------------------------------------------*/
15#ifndef ARCANE_CORE_RANDOM_MERSENNETWISTER_H
16#define ARCANE_CORE_RANDOM_MERSENNETWISTER_H
17/*---------------------------------------------------------------------------*/
18/*---------------------------------------------------------------------------*/
19
20#include "arcane/utils/FatalErrorException.h"
21
22#include "arcane/core/random/RandomGlobal.h"
23
24/*---------------------------------------------------------------------------*/
25/*---------------------------------------------------------------------------*/
26
27namespace Arcane::random
28{
29
30/*---------------------------------------------------------------------------*/
31/*---------------------------------------------------------------------------*/
32
45template<class UIntType, Integer w, Integer n, Integer m, Integer r, UIntType a, Integer u,
46 Integer s, UIntType b, Integer t, UIntType c, Integer l, UIntType val>
48{
49public:
50 typedef UIntType result_type;
51 static const Integer word_size = w;
52 static const Integer state_size = n;
53 static const Integer shift_size = m;
54 static const Integer mask_bits = r;
55 static const UIntType parameter_a = a;
56 static const Integer output_u = u;
57 static const Integer output_s = s;
58 static const UIntType output_b = b;
59 static const Integer output_t = t;
60 static const UIntType output_c = c;
61 static const Integer output_l = l;
62
63 static const bool has_fixed_range = false;
64 /*---------------------------------------------------------------------------*/
65 /*---------------------------------------------------------------------------*/
73 /*---------------------------------------------------------------------------*/
74 /*---------------------------------------------------------------------------*/
81 explicit MersenneTwister(UIntType value)
82 { seed(value); }
83 /*---------------------------------------------------------------------------*/
84 /*---------------------------------------------------------------------------*/
91 template<class It> MersenneTwister(It& first, It last) { seed(first,last); }
92 /*---------------------------------------------------------------------------*/
93 /*---------------------------------------------------------------------------*/
101 template<class Generator>
102 explicit MersenneTwister(Generator & gen) { seed(gen); }
103
104 /*---------------------------------------------------------------------------*/
105 /*---------------------------------------------------------------------------*/
112 void seed() { seed(UIntType(5489)); }
113 /*---------------------------------------------------------------------------*/
114 /*---------------------------------------------------------------------------*/
121 void seed(UIntType value)
122 {
123 // http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/MT2002/emt19937ar.html
124 const UIntType mask = ~0u;
125 x[0] = value & mask;
126 for (i = 1; i < n; i++) {
127 // See Knuth "The Art of Computer Programming" Vol. 2, 3rd ed., page 106
128 x[i] = (1812433253UL * (x[i-1] ^ (x[i-1] >> (w-2))) + i) & mask;
129 }
130 }
131 /*---------------------------------------------------------------------------*/
132 /*---------------------------------------------------------------------------*/
139 void seed(UIntType * state)
140 { for (Integer i=0;i<n;i++) x[i] = state[i];}
141 /*---------------------------------------------------------------------------*/
142 /*---------------------------------------------------------------------------*/
150 template<class Generator>
151 void seed(Generator & gen)
152 {
153 // For GCC, moving this function out-of-line prevents inlining, which may
154 // reduce overall object code size. However, MSVC does not grok
155 // out-of-line definitions of member function templates.
156 for(Integer j = 0; j < n; j++)
157 x[j] = gen();
158 i = n;
159 }
160 /*---------------------------------------------------------------------------*/
161 /*---------------------------------------------------------------------------*/
169 UIntType getState(Integer j)
170 {
171 return x[j];
172 }
173 /*---------------------------------------------------------------------------*/
174 /*---------------------------------------------------------------------------*/
180 result_type min() const { return 0; }
181 /*---------------------------------------------------------------------------*/
182 /*---------------------------------------------------------------------------*/
188 result_type max() const
189 {
190 // avoid "left shift count >= with of type" warning
191 result_type res = 0;
192 for(Integer i = 0; i < w; ++i)
193 res |= (1u << i);
194 return res;
195 }
196 /*---------------------------------------------------------------------------*/
197 /*---------------------------------------------------------------------------*/
198 /* déclaration de l'opérateur () */
199 result_type operator()();
200 /*---------------------------------------------------------------------------*/
201 /*---------------------------------------------------------------------------*/
207 static bool validation(result_type v) { return val == v; }
208 /*---------------------------------------------------------------------------*/
209 /*---------------------------------------------------------------------------*/
215 bool operator==(const MersenneTwister& rhs) const
216 {
217 // Use a member function; Streamable concept not supported.
218 for(Integer j = 0; j < state_size; ++j)
219 if(compute(j) != rhs.compute(j))
220 return false;
221 return true;
222 }
223 /*---------------------------------------------------------------------------*/
224 /*---------------------------------------------------------------------------*/
230 bool operator!=(const MersenneTwister& rhs) const
231 { return !(*this == rhs); }
232
233private:
234 /*---------------------------------------------------------------------------*/
235 /*---------------------------------------------------------------------------*/
242 // returns x(i-n+index), where index is in 0..n-1
243 UIntType compute(UIntType index) const
244 {
245 // equivalent to (i-n+index) % 2n, but doesn't produce negative numbers
246 return x[ (i + n + index) % (2*n) ];
247 }
248 void twist(Integer block);
249 /*---------------------------------------------------------------------------*/
250 /*---------------------------------------------------------------------------*/
251 /* Membres du patron de classe */
252 // state representation: next output is o(x(i))
253 // x[0] ... x[k] x[k+1] ... x[n-1] x[n] ... x[2*n-1] represents
254 // x(i-k) ... x(i) x(i+1) ... x(i-k+n-1) x(i-k-n) ... x[i(i-k-1)]
255 // The goal is to always have x(i-n) ... x(i-1) available for
256 // operator== and save/restore.
257 UIntType x[2*n];
258 Integer i;
259};
260
261/*---------------------------------------------------------------------------*/
262/*---------------------------------------------------------------------------*/
269template<class UIntType, Integer w, Integer n, Integer m, Integer r, UIntType a, Integer u,
270 Integer s, UIntType b, Integer t, UIntType c, Integer l, UIntType val>
272{
273 const UIntType upper_mask = (~0u) << r;
274 const UIntType lower_mask = ~upper_mask;
275
276 if(block == 0) {
277 for(Integer j = n; j < 2*n; j++) {
278 UIntType y = (x[j-n] & upper_mask) | (x[j-(n-1)] & lower_mask);
279 x[j] = x[j-(n-m)] ^ (y >> 1) ^ (y&1 ? a : 0);
280 }
281 } else if (block == 1) {
282 // split loop to avoid costly modulo operations
283 { // extra scope for MSVC brokenness w.r.t. for scope
284 for(Integer j = 0; j < n-m; j++) {
285 UIntType y = (x[j+n] & upper_mask) | (x[j+n+1] & lower_mask);
286 x[j] = x[j+n+m] ^ (y >> 1) ^ (y&1 ? a : 0);
287 }
288 }
289
290 for(Integer j = n-m; j < n-1; j++) {
291 UIntType y = (x[j+n] & upper_mask) | (x[j+n+1] & lower_mask);
292 x[j] = x[j-(n-m)] ^ (y >> 1) ^ (y&1 ? a : 0);
293 }
294 // last iteration
295 UIntType y = (x[2*n-1] & upper_mask) | (x[0] & lower_mask);
296 x[n-1] = x[m-1] ^ (y >> 1) ^ (y&1 ? a : 0);
297 i = 0;
298 }
299}
300
301/*---------------------------------------------------------------------------*/
302/*---------------------------------------------------------------------------*/
309template<class UIntType, Integer w, Integer n, Integer m, Integer r, UIntType a, Integer u,
310 Integer s, UIntType b, Integer t, UIntType c, Integer l, UIntType val>
311inline typename MersenneTwister<UIntType,w,n,m,r,a,u,s,b,t,c,l,val>::result_type
313{
314 if(i == n)
315 twist(0);
316 else if(i >= 2*n)
317 twist(1);
318 // Step 4
319 UIntType z = x[i];
320 ++i;
321 z ^= (z >> u);
322 z ^= ((z << s) & b);
323 z ^= ((z << t) & c);
324 z ^= (z >> l);
325 return z;
326}
327
328/*---------------------------------------------------------------------------*/
329/*---------------------------------------------------------------------------*/
330/* Définition de la classe mt11213b*/
331typedef MersenneTwister<UInt32,32,351,175,19,0xccab8ee7,11,
332 7,0x31b6ab00,15,0xffe50000,17, 0xa37d3c92> Mt11213b;
333/*---------------------------------------------------------------------------*/
334/*---------------------------------------------------------------------------*/
335/* Définition de la classe Mt19937*/
336typedef MersenneTwister<UInt32,32,624,397,31,0x9908b0df,11,
337 7,0x9d2c5680,15,0xefc60000,18, 3346425566U> Mt19937;
338
339}
340
341/*---------------------------------------------------------------------------*/
342/*---------------------------------------------------------------------------*/
343
344#endif // ARCANE_RANDOM_MERSENNE_TWISTER_H
Patron de classe MersenneTwister.
static bool validation(result_type v)
Fonction de validation (je ne sais pas trop a quoi elle sert!)
MersenneTwister()
Constructeur avec initialisation de la graine à partir de la méthode seed()
MersenneTwister(UIntType value)
Constructeur avec initialisation du tableau de graines à partir de la valeur value....
bool operator!=(const MersenneTwister &rhs) const
Surdéfinition de l'opérateur !=.
void seed(UIntType *state)
Initialisation du tableau de graines à partir du tableau state. state doit être un tableau de n éléme...
void seed(UIntType value)
Initialisation du tableau de graines à partir de la valeur value. Le tableau de graines de ce générat...
MersenneTwister(Generator &gen)
Constructeur avec initialisation du tableau de graines à partir du générateur gen....
result_type max() const
max() retourne la valeur maximum possible d'une séquence.
void seed(Generator &gen)
Initialisation du tableau de graines à partir du générateur gen. gen est une classe qui doit contenir...
UIntType getState(Integer j)
Méthode qui retourne l'état du générateur pour l'index j. L'état complet du générateur est donnée par...
bool operator==(const MersenneTwister &rhs) const
Surdéfinition de l'opérateur ==.
void twist(Integer block)
Réalisation de l'opération "twist" associée au Mersenne Twister. L'état du générateur est modifié.
UIntType compute(UIntType index) const
Méthode privée qui retourne l'état du générateur pour l'index index.
MersenneTwister(It &first, It last)
Constructeur avec initialisation du tableau de graines à partir de la méthode seed(first,...
void seed()
Initialisation du tableau de graines. L'appel à la méthode seed(5489) est réalisé.
result_type operator()()
Surdéfinition de l'opérateur () qui retourne la valeur pseudo aléatoire du générateur....
result_type min() const
min() retourne la valeur minimum possible d'une séquence.
std::uint32_t UInt32
Type entier non signé sur 32 bits.
Int32 Integer
Type représentant un entier.