cachedRandom.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 2011-2015 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
26 #include "cachedRandom.H"
27 #include "OSspecific.H"
28 #include "PstreamReduceOps.H"
29 
30 // * * * * * * * * * * * * * private Member Functions * * * * * * * * * * * //
31 
32 Foam::scalar Foam::cachedRandom::scalar01()
33 {
34  if (sampleI_ < 0)
35  {
36  return osRandomDouble();
37  }
38 
39  if (sampleI_ == samples_.size() - 1)
40  {
41  scalar s = samples_[sampleI_];
42  sampleI_ = 0;
43  return s;
44  }
45  else
46  {
47  scalar s = samples_[sampleI_];
48  sampleI_++;
49  return s;
50  }
51 }
52 
53 
54 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
55 
57 :
58  seed_(1),
59  samples_(0),
60  sampleI_(-1),
61  hasGaussSample_(false),
62  gaussSample_(0)
63 {
64  if (seed > 1)
65  {
66  seed_ = seed;
67  }
68 
69  // Samples will be cached if count > 0
70  if (count > 0)
71  {
72  samples_.setSize(count);
73  sampleI_ = 0;
74  }
75 
76  // Initialise samples
77  osRandomSeed(seed_);
78  forAll(samples_, i)
79  {
80  samples_[i] = osRandomDouble();
81  }
82 }
83 
84 
85 Foam::cachedRandom::cachedRandom(const cachedRandom& cr, const bool reset)
86 :
87  seed_(cr.seed_),
88  samples_(cr.samples_),
89  sampleI_(cr.sampleI_),
90  hasGaussSample_(cr.hasGaussSample_),
91  gaussSample_(cr.gaussSample_)
92 {
93  if (reset)
94  {
95  hasGaussSample_ = false;
96  gaussSample_ = 0;
97  }
98  if (sampleI_ == -1)
99  {
100  WarningIn
101  (
102  "Foam::cachedRandom::cachedRandom(const cachedRandom& cr)"
103  ) << "Copy constructor called, but samples not being cached. "
104  << "This may lead to non-repeatable behaviour" << endl;
105 
106  osRandomSeed(seed_);
107  }
108 
109  if (reset && samples_.size())
110  {
111  sampleI_ = 0;
112  }
113 }
114 
115 
116 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
117 
119 {}
120 
121 
122 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
123 
124 template<>
126 {
127  return scalar01();
128 }
129 
130 
131 template<>
133 {
134  return round(scalar01());
135 }
136 
137 
138 template<>
140 {
141  if (hasGaussSample_)
142  {
143  hasGaussSample_ = false;
144  return gaussSample_;
145  }
146  else
147  {
148  scalar rsq, v1, v2;
149  do
150  {
151  v1 = 2*scalar01() - 1;
152  v2 = 2*scalar01() - 1;
153  rsq = sqr(v1) + sqr(v2);
154  } while (rsq >= 1 || rsq == 0);
155 
156  scalar fac = sqrt(-2*log(rsq)/rsq);
157 
158  gaussSample_ = v1*fac;
159  hasGaussSample_ = true;
160 
161  return v2*fac;
162  }
163 }
164 
165 
166 template<>
168 {
169  return round(GaussNormal<scalar>());
170 }
171 
172 
173 template<>
174 Foam::scalar Foam::cachedRandom::position
175 (
176  const scalar& start,
177  const scalar& end
178 )
179 {
180  return start + scalar01()*(end - start);
181 }
182 
183 
184 template<>
186 {
187  return start + round(scalar01()*(end - start));
188 }
189 
190 
191 template<>
193 {
194  scalar value = -GREAT;
195 
196  if (Pstream::master())
197  {
198  value = scalar01();
199  }
200 
201  Pstream::scatter(value);
202 
203  return value;
204 }
205 
206 
207 template<>
209 {
210  scalar value = -GREAT;
211 
212  if (Pstream::master())
213  {
214  value = scalar01();
215  }
216 
217  Pstream::scatter(value);
218 
219  return round(value);
220 }
221 
222 
223 template<>
225 {
226  scalar value = -GREAT;
227 
228  if (Pstream::master())
229  {
230  value = GaussNormal<scalar>();
231  }
232 
233  Pstream::scatter(value);
234 
235  return value;
236 }
237 
238 
239 template<>
241 {
242  scalar value = -GREAT;
243 
244  if (Pstream::master())
245  {
246  value = GaussNormal<scalar>();
247  }
248 
249  Pstream::scatter(value);
250 
251  return round(value);
252 }
253 
254 
255 template<>
257 (
258  const scalar& start,
259  const scalar& end
260 )
261 {
262  scalar value = -GREAT;
263 
264  if (Pstream::master())
265  {
266  value = scalar01()*(end - start);
267  }
268 
269  Pstream::scatter(value);
270 
271  return start + value;
272 }
273 
274 
275 template<>
277 (
278  const label& start,
279  const label& end
280 )
281 {
282  label value = labelMin;
283 
284  if (Pstream::master())
285  {
286  value = round(scalar01()*(end - start));
287  }
288 
289  Pstream::scatter(value);
290 
291  return start + value;
292 }
293 
294 
295 // ************************************************************************* //
dimensionedScalar sqrt(const dimensionedScalar &ds)
Type globalSample01()
Return a sample whose components lie in the range 0-1.
gmvFile<< "tracers "<< particles.size()<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().x()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().y()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().z()<< " ";}gmvFile<< nl;forAll(lagrangianScalarNames, i){word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject( name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE ))
scalar osRandomDouble()
Return random double precision (uniform distribution between 0 and 1)
Definition: POSIX.C:1311
Functions used by OpenFOAM that are specific to POSIX compliant operating systems and need to be repl...
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
Type globalPosition(const Type &start, const Type &end)
Return a sample between start and end.
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:76
dimensionedScalar log(const dimensionedScalar &ds)
void setSize(const label)
Reset size of List.
Definition: List.C:318
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
#define WarningIn(functionName)
Report a warning using Foam::Warning.
#define forAll(list, i)
Definition: UList.H:421
static void scatter(const List< commsStruct > &comms, T &Value, const int tag, const label comm)
Scatter data. Distribute without modification. Reverse of gather.
static const label labelMin
Definition: label.H:61
void osRandomSeed(const label seed)
Seed random number generator.
Definition: POSIX.C:1291
Random number generator.
Definition: cachedRandom.H:63
Type sample01()
Return a sample whose components lie in the range 0-1.
static bool master(const label communicator=0)
Am I the master process.
Definition: UPstream.H:398
label seed() const
Return const access to the initial random number seed.
Definition: cachedRandomI.H:30
Type position(const Type &start, const Type &end)
Return a sample between start and end.
Type globalGaussNormal()
Return a sample whose components are normally distributed.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
Type GaussNormal()
Return a sample whose components are normally distributed.
cachedRandom(const label seed, const label count)
Construct given seed and sample count.
Definition: cachedRandom.C:56