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