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