CASToR  3.2
Tomographic Reconstruction (PET/SPECT/CT)
code/src/management/sRandomNumberGenerator.cc
Go to the documentation of this file.
1 
10 #include "sRandomNumberGenerator.hh"
11 
13 
14 // =====================================================================
15 // ---------------------------------------------------------------------
16 // ---------------------------------------------------------------------
17 // =====================================================================
22 {
23  mp_Instance = NULL;
24  m_verbose = 5;
25  mp_Engines.clear();
26  mp_extraEngines.clear();
27  m_seed = -1;
28 }
29 
30 
31 
32 // =====================================================================
33 // ---------------------------------------------------------------------
34 // ---------------------------------------------------------------------
35 // =====================================================================
36 /*
37  \brief Destructor of sRandomNumberGenerator. Do nothing by default
38 */
40 {
41 }
42 
43 
44 
45 // =====================================================================
46 // ---------------------------------------------------------------------
47 // ---------------------------------------------------------------------
48 // =====================================================================
49 
50 int sRandomNumberGenerator::Initialize(int a_nbThreads, int a_nbExtra)
51 {
52  if(m_verbose >=3) Cout("sRandomNumberGenerator::Initialize ..."<< endl);
53 
54  int64_t seed = -1;
55 
56  // call a real random generator, only once because this operating system
57  // random source might be limited in size and time
58  std::random_device rd;
59 
60  // save the initial seed
61  m_seed = rd();
62 
63  #ifdef CASTOR_MPI
64  int mpi_rank_temp = 0;
65  int mpi_size_temp = 0;
66  MPI_Comm_rank(MPI_COMM_WORLD, &mpi_rank_temp);
67  MPI_Comm_size(MPI_COMM_WORLD, &mpi_size_temp);
68 
69  // helper variable
70  int64_t temp_seed = -1;
71 
72  // if more than one MPI instance, generate seeds for each instance
73  // and dispatch them, otherwise keep the initial seed
74  if (mpi_size_temp>1)
75  {
76  if (mpi_rank_temp==0)
77  {
78  Engine mpi_generator(m_seed);
79  for (int p=0; p<mpi_size_temp; p++)
80  {
81  m_seed = mpi_generator();
82 
83  Cout("sRandomNumberGenerator::Seed for rank " << p << " is " << m_seed << endl);
84 
85  if (p==0)
86  temp_seed = m_seed;
87  else
88  MPI_Send(&m_seed, 1, MPI_LONG, p, 0, MPI_COMM_WORLD);
89  }
90  m_seed = temp_seed;
91  }
92  else
93  {
94  MPI_Status* status = NULL;
95  MPI_Recv(&m_seed, 1, MPI_LONG, 0, 0, MPI_COMM_WORLD, status);
96  }
97  // wait for all the processes to have their seeds
98  MPI_Barrier(MPI_COMM_WORLD);
99  }
100  #endif
101 
102  // seed the initial temporary PRNG, which will generate seeds for actual PRNGs (threaded and non-threaded)
103  Engine initialGenerator(m_seed);
104 
105  for (int th=0; th<a_nbThreads; th++)
106  {
107  seed = initialGenerator();
108  if (m_verbose>=3) Cout("sRandomNumberGenerator::Initialize() -> Seed for thread " << th << " : " << seed << endl);
109  mp_Engines.push_back(Engine(seed));
110  }
111 
112  for(int ex=0; ex<a_nbExtra; ex++)
113  {
114  seed = initialGenerator();
115  if(m_verbose>=3) Cout("sRandomNumberGenerator::Initialize() - >Seed for additional nonthreaded generator "<<ex<<" : " << seed << endl);
116  mp_extraEngines.push_back(Engine(seed));
117  }
118 
119  return 0;
120 }
121 
122 
123 
124 // =====================================================================
125 // ---------------------------------------------------------------------
126 // ---------------------------------------------------------------------
127 // =====================================================================
128 
129 int sRandomNumberGenerator::Initialize(int64_t a_seed, int a_nbThreads, int a_nbExtra)
130 {
131  if(m_verbose >=3) Cout("sRandomNumberGenerator::Initialize with provided seed "<<a_seed<< endl);
132 
133  // temporary variable
134  int64_t seed = -1;
135 
136  if (a_seed<0)
137  {
138  Cout("***** sRandomNumberGenerator::Initialize()-> Error : seed for RNG should be >=0 !" << endl);
139  return 1;
140  }
141 
142  // save the initial seed
143  m_seed = a_seed;
144 
145  #ifdef CASTOR_MPI
146  int mpi_rank_temp = 0;
147  int mpi_size_temp = 0;
148  MPI_Comm_rank(MPI_COMM_WORLD, &mpi_rank_temp);
149  MPI_Comm_size(MPI_COMM_WORLD, &mpi_size_temp);
150 
151  // helper variable
152  int64_t temp_seed = -1;
153 
154  // if more than one MPI instance, generate seeds for each instance
155  // and dispatch them, otherwise keep the initial seed
156  if (mpi_size_temp>1)
157  {
158  if (mpi_rank_temp==0)
159  {
160  Engine mpi_generator(a_seed);
161  for (int p=0; p<mpi_size_temp; p++)
162  {
163  m_seed = mpi_generator();
164 
165  Cout("sRandomNumberGenerator::Seed for rank "<<p<<" is "<<m_seed<<""<< endl);
166 
167  if (p==0)
168  temp_seed = m_seed;
169  else
170  MPI_Send(&m_seed, 1, MPI_LONG, p, 0, MPI_COMM_WORLD);
171  }
172  m_seed = temp_seed;
173  }
174  else
175  {
176  MPI_Status* status = NULL;
177  MPI_Recv(&m_seed, 1, MPI_LONG, 0, 0, MPI_COMM_WORLD, status);
178  }
179  MPI_Barrier(MPI_COMM_WORLD);
180  }
181  #endif
182 
183  // seed the initial temporary PRNG, which will generate seeds for actual PRNGs (threaded and non-threaded)
184  Engine initialGenerator(m_seed);
185 
186  for(int th=0; th<a_nbThreads; th++)
187  {
188  seed = initialGenerator();
189  mp_Engines.push_back(Engine(seed));
190  if (m_verbose>=3) Cout("sRandomNumberGenerator::Initialize() -> Seed for thread " << th << " : " << seed << endl);
191  }
192 
193  for(int ex=0; ex<a_nbExtra; ex++)
194  {
195  seed = initialGenerator();
196  mp_extraEngines.push_back(Engine(seed));
197  if (m_verbose>=3) Cout("sRandomNumberGenerator::Initialize() -> Seed for additional nonthreaded random generator "<<ex<< " : " << seed << endl);
198  }
199 
200  return 0;
201 }
202 
203 
204 
205 // =====================================================================
206 // ---------------------------------------------------------------------
207 // ---------------------------------------------------------------------
208 // =====================================================================
209 
211 {
212  #ifdef CASTOR_VERBOSE
213  if (m_verbose >=4) Cout("sRandomNumberGenerator::GenerateRdmNber..." << endl);
214  #endif
215 
216  int id=0;
217  #ifdef CASTOR_OMP
218  id = omp_get_thread_num();
219  #endif
220  return mp_Distribution(mp_Engines[id]);
221 }
222 
223 
224 
225 // =====================================================================
226 // ---------------------------------------------------------------------
227 // ---------------------------------------------------------------------
228 // =====================================================================
229 
231 {
232  #ifdef CASTOR_VERBOSE
233  if(m_verbose >=4) Cout("sRandomNumberGenerator::GenerateExtraRdmNber..."<< endl);
234  #endif
235 
236  return mp_Distribution(mp_extraEngines[a_nb]);
237 }
238 
239 
240 // =====================================================================
241 // ---------------------------------------------------------------------
242 // ---------------------------------------------------------------------
243 // =====================================================================
244 
246 {
247  #ifdef CASTOR_VERBOSE
248  if(m_verbose >=4) Cout("sRandomNumberGenerator::GetExtraGenerator..."<< endl);
249  #endif
250 
251  return mp_extraEngines[a_nb];
252 }
253 
254 
255 // =====================================================================
256 // ---------------------------------------------------------------------
257 // ---------------------------------------------------------------------
258 // =====================================================================
259 
261 {
262  #ifdef CASTOR_VERBOSE
263  if(m_verbose >=4) Cout("sRandomNumberGenerator::GetGenerator..."<< endl);
264  #endif
265 
266  int id=0;
267  #ifdef CASTOR_OMP
268  id = omp_get_thread_num();
269  #endif
270 
271  return mp_Engines[id];
272 }
Engine & GetExtraGenerator(int a_nb=0)
Get the not thread safe additional random generator, for use various sequential parts of an otherwise...
HPFLTNB GenerateExtraRdmNber(int a_nb=0)
Generate a random number using the specified additional not thread safe random generator, for use in sequential parts of an otherwise multithreaded code.
Engine & GetGenerator()
Get the random generator for the current thread.
int Initialize(int a_nbThreads, int a_nbExtra)
Instantiate pseudo random number generators, one per thread by default, and additional extra ones if ...
sRandomNumberGenerator()
Constructor of sRandomNumberGenerator. Do nothing by default as it is a singleton clasee...
HPFLTNB GenerateRdmNber()
Generate a random number for the thread which index is recovered from the OpenMP function.
Singleton class that generate a thread-safe random generator number for openMP As singleton...
~sRandomNumberGenerator()
Destructor of sRandomNumberGenerator. Do nothing by default.
#define Cout(MESSAGE)