1
2 /*
3 This file contains routines for interfacing to random number generators.
4 This provides more than just an interface to some system random number
5 generator:
6
7 Numbers can be shuffled for use as random tuples
8
9 Multiple random number generators may be used
10
11 We are still not sure what interface we want here. There should be
12 one to reinitialize and set the seed.
13 */
14
15 #include <../src/sys/classes/random/randomimpl.h> /*I "petscsys.h" I*/
16 #include <petscviewer.h>
17
18 /* Logging support */
19 PetscClassId PETSC_RANDOM_CLASSID;
20
21 /*@
22 PetscRandomDestroy - Destroys a context that has been formed by
23 PetscRandomCreate().
24
25 Collective on PetscRandom
26
27 Intput Parameter:
28 . r - the random number generator context
29
30 Level: intermediate
31
32 .seealso: PetscRandomGetValue(), PetscRandomCreate(), VecSetRandom()
33 @*/
PetscRandomDestroy(PetscRandom * r)34 PetscErrorCode PetscRandomDestroy(PetscRandom *r)
35 {
36 PetscErrorCode ierr;
37
38 PetscFunctionBegin;
39 if (!*r) PetscFunctionReturn(0);
40 PetscValidHeaderSpecific(*r,PETSC_RANDOM_CLASSID,1);
41 if (--((PetscObject)(*r))->refct > 0) {*r = NULL; PetscFunctionReturn(0);}
42 if ((*r)->ops->destroy) {
43 ierr = (*(*r)->ops->destroy)(*r);CHKERRQ(ierr);
44 }
45 ierr = PetscHeaderDestroy(r);CHKERRQ(ierr);
46 PetscFunctionReturn(0);
47 }
48
49
50 /*@C
51 PetscRandomGetSeed - Gets the random seed.
52
53 Not collective
54
55 Input Parameters:
56 . r - The random number generator context
57
58 Output Parameter:
59 . seed - The random seed
60
61 Level: intermediate
62
63 .seealso: PetscRandomCreate(), PetscRandomSetSeed(), PetscRandomSeed()
64 @*/
PetscRandomGetSeed(PetscRandom r,unsigned long * seed)65 PetscErrorCode PetscRandomGetSeed(PetscRandom r,unsigned long *seed)
66 {
67 PetscFunctionBegin;
68 PetscValidHeaderSpecific(r,PETSC_RANDOM_CLASSID,1);
69 if (seed) {
70 PetscValidPointer(seed,2);
71 *seed = r->seed;
72 }
73 PetscFunctionReturn(0);
74 }
75
76 /*@C
77 PetscRandomSetSeed - Sets the random seed. You MUST call PetscRandomSeed() after this call to have the new seed used.
78
79 Not collective
80
81 Input Parameters:
82 + r - The random number generator context
83 - seed - The random seed
84
85 Level: intermediate
86
87 Usage:
88 PetscRandomSetSeed(r,a positive integer);
89 PetscRandomSeed(r); PetscRandomGetValue() will now start with the new seed.
90
91 PetscRandomSeed(r) without a call to PetscRandomSetSeed() re-initializes
92 the seed. The random numbers generated will be the same as before.
93
94 .seealso: PetscRandomCreate(), PetscRandomGetSeed(), PetscRandomSeed()
95 @*/
PetscRandomSetSeed(PetscRandom r,unsigned long seed)96 PetscErrorCode PetscRandomSetSeed(PetscRandom r,unsigned long seed)
97 {
98 PetscErrorCode ierr;
99
100 PetscFunctionBegin;
101 PetscValidHeaderSpecific(r,PETSC_RANDOM_CLASSID,1);
102 r->seed = seed;
103 ierr = PetscInfo1(NULL,"Setting seed to %d\n",(int)seed);CHKERRQ(ierr);
104 PetscFunctionReturn(0);
105 }
106
107 /* ------------------------------------------------------------------- */
108 /*
109 PetscRandomSetTypeFromOptions_Private - Sets the type of random generator from user options. Defaults to type PETSCRAND48 or PETSCRAND.
110
111 Collective on PetscRandom
112
113 Input Parameter:
114 . rnd - The random number generator context
115
116 Level: intermediate
117
118 .seealso: PetscRandomSetFromOptions(), PetscRandomSetType()
119 */
PetscRandomSetTypeFromOptions_Private(PetscOptionItems * PetscOptionsObject,PetscRandom rnd)120 static PetscErrorCode PetscRandomSetTypeFromOptions_Private(PetscOptionItems *PetscOptionsObject,PetscRandom rnd)
121 {
122 PetscBool opt;
123 const char *defaultType;
124 char typeName[256];
125 PetscErrorCode ierr;
126
127 PetscFunctionBegin;
128 if (((PetscObject)rnd)->type_name) {
129 defaultType = ((PetscObject)rnd)->type_name;
130 } else {
131 defaultType = PETSCRANDER48;
132 }
133
134 ierr = PetscRandomRegisterAll();CHKERRQ(ierr);
135 ierr = PetscOptionsFList("-random_type","PetscRandom type","PetscRandomSetType",PetscRandomList,defaultType,typeName,256,&opt);CHKERRQ(ierr);
136 if (opt) {
137 ierr = PetscRandomSetType(rnd, typeName);CHKERRQ(ierr);
138 } else {
139 ierr = PetscRandomSetType(rnd, defaultType);CHKERRQ(ierr);
140 }
141 PetscFunctionReturn(0);
142 }
143
144 /*@
145 PetscRandomSetFromOptions - Configures the random number generator from the options database.
146
147 Collective on PetscRandom
148
149 Input Parameter:
150 . rnd - The random number generator context
151
152 Options Database:
153 + -random_seed <integer> - provide a seed to the random number generater
154 - -random_no_imaginary_part - makes the imaginary part of the random number zero, this is useful when you want the
155 same code to produce the same result when run with real numbers or complex numbers for regression testing purposes
156
157 Notes:
158 To see all options, run your program with the -help option.
159 Must be called after PetscRandomCreate() but before the rnd is used.
160
161 Level: beginner
162
163 .seealso: PetscRandomCreate(), PetscRandomSetType()
164 @*/
PetscRandomSetFromOptions(PetscRandom rnd)165 PetscErrorCode PetscRandomSetFromOptions(PetscRandom rnd)
166 {
167 PetscErrorCode ierr;
168 PetscBool set,noimaginary = PETSC_FALSE;
169 PetscInt seed;
170
171 PetscFunctionBegin;
172 PetscValidHeaderSpecific(rnd,PETSC_RANDOM_CLASSID,1);
173
174 ierr = PetscObjectOptionsBegin((PetscObject)rnd);CHKERRQ(ierr);
175
176 /* Handle PetscRandom type options */
177 ierr = PetscRandomSetTypeFromOptions_Private(PetscOptionsObject,rnd);CHKERRQ(ierr);
178
179 /* Handle specific random generator's options */
180 if (rnd->ops->setfromoptions) {
181 ierr = (*rnd->ops->setfromoptions)(PetscOptionsObject,rnd);CHKERRQ(ierr);
182 }
183 ierr = PetscOptionsInt("-random_seed","Seed to use to generate random numbers","PetscRandomSetSeed",0,&seed,&set);CHKERRQ(ierr);
184 if (set) {
185 ierr = PetscRandomSetSeed(rnd,(unsigned long int)seed);CHKERRQ(ierr);
186 ierr = PetscRandomSeed(rnd);CHKERRQ(ierr);
187 }
188 ierr = PetscOptionsBool("-random_no_imaginary_part","The imaginary part of the random number will be zero","PetscRandomSetInterval",noimaginary,&noimaginary,&set);CHKERRQ(ierr);
189 #if defined(PETSC_HAVE_COMPLEX)
190 if (set) {
191 if (noimaginary) {
192 PetscScalar low,high;
193 ierr = PetscRandomGetInterval(rnd,&low,&high);CHKERRQ(ierr);
194 low = low - PetscImaginaryPart(low);
195 high = high - PetscImaginaryPart(high);
196 ierr = PetscRandomSetInterval(rnd,low,high);CHKERRQ(ierr);
197 }
198 }
199 #endif
200 ierr = PetscOptionsEnd();CHKERRQ(ierr);
201 ierr = PetscRandomViewFromOptions(rnd,NULL, "-random_view");CHKERRQ(ierr);
202 PetscFunctionReturn(0);
203 }
204
205 #if defined(PETSC_HAVE_SAWS)
206 #include <petscviewersaws.h>
207 #endif
208
209 /*@C
210 PetscRandomViewFromOptions - View from Options
211
212 Collective on PetscRandom
213
214 Input Parameters:
215 + A - the random number generator context
216 . obj - Optional object
217 - name - command line option
218
219 Level: intermediate
220 .seealso: PetscRandom, PetscRandomView, PetscObjectViewFromOptions(), PetscRandomCreate()
221 @*/
PetscRandomViewFromOptions(PetscRandom A,PetscObject obj,const char name[])222 PetscErrorCode PetscRandomViewFromOptions(PetscRandom A,PetscObject obj,const char name[])
223 {
224 PetscErrorCode ierr;
225
226 PetscFunctionBegin;
227 PetscValidHeaderSpecific(A,PETSC_RANDOM_CLASSID,1);
228 ierr = PetscObjectViewFromOptions((PetscObject)A,obj,name);CHKERRQ(ierr);
229 PetscFunctionReturn(0);
230 }
231
232 /*@C
233 PetscRandomView - Views a random number generator object.
234
235 Collective on PetscRandom
236
237 Input Parameters:
238 + rnd - The random number generator context
239 - viewer - an optional visualization context
240
241 Notes:
242 The available visualization contexts include
243 + PETSC_VIEWER_STDOUT_SELF - standard output (default)
244 - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
245 output where only the first processor opens
246 the file. All other processors send their
247 data to the first processor to print.
248
249 You can change the format the vector is printed using the
250 option PetscViewerPushFormat().
251
252 Level: beginner
253
254 .seealso: PetscRealView(), PetscScalarView(), PetscIntView()
255 @*/
PetscRandomView(PetscRandom rnd,PetscViewer viewer)256 PetscErrorCode PetscRandomView(PetscRandom rnd,PetscViewer viewer)
257 {
258 PetscErrorCode ierr;
259 PetscBool iascii;
260 #if defined(PETSC_HAVE_SAWS)
261 PetscBool issaws;
262 #endif
263
264 PetscFunctionBegin;
265 PetscValidHeaderSpecific(rnd,PETSC_RANDOM_CLASSID,1);
266 PetscValidType(rnd,1);
267 if (!viewer) {
268 ierr = PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)rnd),&viewer);CHKERRQ(ierr);
269 }
270 PetscValidHeaderSpecific(viewer,PETSC_VIEWER_CLASSID,2);
271 PetscCheckSameComm(rnd,1,viewer,2);
272 ierr = PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);CHKERRQ(ierr);
273 #if defined(PETSC_HAVE_SAWS)
274 ierr = PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERSAWS,&issaws);CHKERRQ(ierr);
275 #endif
276 if (iascii) {
277 PetscMPIInt rank;
278 ierr = PetscObjectPrintClassNamePrefixType((PetscObject)rnd,viewer);CHKERRQ(ierr);
279 ierr = MPI_Comm_rank(PetscObjectComm((PetscObject)rnd),&rank);CHKERRQ(ierr);
280 ierr = PetscViewerASCIIPushSynchronized(viewer);CHKERRQ(ierr);
281 ierr = PetscViewerASCIISynchronizedPrintf(viewer,"[%d] Random type %s, seed %lu\n",rank,((PetscObject)rnd)->type_name,rnd->seed);CHKERRQ(ierr);
282 ierr = PetscViewerFlush(viewer);CHKERRQ(ierr);
283 ierr = PetscViewerASCIIPopSynchronized(viewer);CHKERRQ(ierr);
284 #if defined(PETSC_HAVE_SAWS)
285 } else if (issaws) {
286 PetscMPIInt rank;
287 const char *name;
288
289 ierr = PetscObjectGetName((PetscObject)rnd,&name);CHKERRQ(ierr);
290 ierr = MPI_Comm_rank(PETSC_COMM_WORLD,&rank);CHKERRQ(ierr);
291 if (!((PetscObject)rnd)->amsmem && !rank) {
292 char dir[1024];
293
294 ierr = PetscObjectViewSAWs((PetscObject)rnd,viewer);CHKERRQ(ierr);
295 ierr = PetscSNPrintf(dir,1024,"/PETSc/Objects/%s/Low",name);CHKERRQ(ierr);
296 PetscStackCallSAWs(SAWs_Register,(dir,&rnd->low,1,SAWs_READ,SAWs_DOUBLE));
297 }
298 #endif
299 }
300 PetscFunctionReturn(0);
301 }
302
303 /*@
304 PetscRandomCreate - Creates a context for generating random numbers,
305 and initializes the random-number generator.
306
307 Collective
308
309 Input Parameters:
310 . comm - MPI communicator
311
312 Output Parameter:
313 . r - the random number generator context
314
315 Level: intermediate
316
317 Notes:
318 The random type has to be set by PetscRandomSetType().
319
320 This is only a primative "parallel" random number generator, it should NOT
321 be used for sophisticated parallel Monte Carlo methods since it will very likely
322 not have the correct statistics across processors. You can provide your own
323 parallel generator using PetscRandomRegister();
324
325 If you create a PetscRandom() using PETSC_COMM_SELF on several processors then
326 the SAME random numbers will be generated on all those processors. Use PETSC_COMM_WORLD
327 or the appropriate parallel communicator to eliminate this issue.
328
329 Use VecSetRandom() to set the elements of a vector to random numbers.
330
331 Example of Usage:
332 .vb
333 PetscRandomCreate(PETSC_COMM_SELF,&r);
334 PetscRandomSetType(r,PETSCRAND48);
335 PetscRandomGetValue(r,&value1);
336 PetscRandomGetValueReal(r,&value2);
337 PetscRandomDestroy(&r);
338 .ve
339
340 .seealso: PetscRandomSetType(), PetscRandomGetValue(), PetscRandomGetValueReal(), PetscRandomSetInterval(),
341 PetscRandomDestroy(), VecSetRandom(), PetscRandomType
342 @*/
343
PetscRandomCreate(MPI_Comm comm,PetscRandom * r)344 PetscErrorCode PetscRandomCreate(MPI_Comm comm,PetscRandom *r)
345 {
346 PetscRandom rr;
347 PetscErrorCode ierr;
348 PetscMPIInt rank;
349
350 PetscFunctionBegin;
351 PetscValidPointer(r,3);
352 *r = NULL;
353 ierr = PetscRandomInitializePackage();CHKERRQ(ierr);
354
355 ierr = PetscHeaderCreate(rr,PETSC_RANDOM_CLASSID,"PetscRandom","Random number generator","Sys",comm,PetscRandomDestroy,PetscRandomView);CHKERRQ(ierr);
356
357 ierr = MPI_Comm_rank(comm,&rank);CHKERRQ(ierr);
358
359 rr->data = NULL;
360 rr->low = 0.0;
361 rr->width = 1.0;
362 rr->iset = PETSC_FALSE;
363 rr->seed = 0x12345678 + 76543*rank;
364 ierr = PetscRandomSetType(rr,PETSCRANDER48);CHKERRQ(ierr);
365 *r = rr;
366 PetscFunctionReturn(0);
367 }
368
369 /*@
370 PetscRandomSeed - Seed the generator.
371
372 Not collective
373
374 Input Parameters:
375 . r - The random number generator context
376
377 Level: intermediate
378
379 Usage:
380 PetscRandomSetSeed(r,a positive integer);
381 PetscRandomSeed(r); PetscRandomGetValue() will now start with the new seed.
382
383 PetscRandomSeed(r) without a call to PetscRandomSetSeed() re-initializes
384 the seed. The random numbers generated will be the same as before.
385
386 .seealso: PetscRandomCreate(), PetscRandomGetSeed(), PetscRandomSetSeed()
387 @*/
PetscRandomSeed(PetscRandom r)388 PetscErrorCode PetscRandomSeed(PetscRandom r)
389 {
390 PetscErrorCode ierr;
391
392 PetscFunctionBegin;
393 PetscValidHeaderSpecific(r,PETSC_RANDOM_CLASSID,1);
394 PetscValidType(r,1);
395
396 ierr = (*r->ops->seed)(r);CHKERRQ(ierr);
397 ierr = PetscObjectStateIncrease((PetscObject)r);CHKERRQ(ierr);
398 PetscFunctionReturn(0);
399 }
400
401