xref: /petsc/src/sys/classes/random/interface/randomc.c (revision af0996ce37bc06907c37d8d91773840993d61e62)
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(PetscOptions *PetscOptionsObject,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   ierr = PetscRandomRegisterAll();CHKERRQ(ierr);
149   ierr = PetscOptionsFList("-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(PetscOptionsObject,rnd);CHKERRQ(ierr);
192 
193   /* Handle specific random generator's options */
194   if (rnd->ops->setfromoptions) {
195     ierr = (*rnd->ops->setfromoptions)(PetscOptionsObject,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      issaws;
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,&issaws);CHKERRQ(ierr);
255 #endif
256   if (iascii) {
257     PetscMPIInt rank;
258     ierr = PetscObjectPrintClassNamePrefixType((PetscObject)rnd,viewer);CHKERRQ(ierr);
259     ierr = MPI_Comm_rank(PetscObjectComm((PetscObject)rnd),&rank);CHKERRQ(ierr);
260     ierr = PetscViewerASCIISynchronizedAllow(viewer,PETSC_TRUE);CHKERRQ(ierr);
261     ierr = PetscViewerASCIISynchronizedPrintf(viewer,"[%D] Random type %s, seed %D\n",rank,((PetscObject)rnd)->type_name,rnd->seed);CHKERRQ(ierr);
262     ierr = PetscViewerFlush(viewer);CHKERRQ(ierr);
263     ierr = PetscViewerASCIISynchronizedAllow(viewer,PETSC_FALSE);CHKERRQ(ierr);
264 #if defined(PETSC_HAVE_SAWS)
265   } else if (issaws) {
266     PetscMPIInt rank;
267     const char  *name;
268 
269     ierr = PetscObjectGetName((PetscObject)rnd,&name);CHKERRQ(ierr);
270     ierr = MPI_Comm_rank(PETSC_COMM_WORLD,&rank);CHKERRQ(ierr);
271     if (!((PetscObject)rnd)->amsmem && !rank) {
272       char       dir[1024];
273 
274       ierr = PetscObjectViewSAWs((PetscObject)rnd,viewer);CHKERRQ(ierr);
275       ierr = PetscSNPrintf(dir,1024,"/PETSc/Objects/%s/Low",name);CHKERRQ(ierr);
276       PetscStackCallSAWs(SAWs_Register,(dir,&rnd->low,1,SAWs_READ,SAWs_DOUBLE));
277     }
278 #endif
279   }
280   PetscFunctionReturn(0);
281 }
282 
283 #undef __FUNCT__
284 #define __FUNCT__ "PetscRandomCreate"
285 /*@
286    PetscRandomCreate - Creates a context for generating random numbers,
287    and initializes the random-number generator.
288 
289    Collective on MPI_Comm
290 
291    Input Parameters:
292 +  comm - MPI communicator
293 
294    Output Parameter:
295 .  r  - the random number generator context
296 
297    Level: intermediate
298 
299    Notes:
300    The random type has to be set by PetscRandomSetType().
301 
302    This is only a primative "parallel" random number generator, it should NOT
303    be used for sophisticated parallel Monte Carlo methods since it will very likely
304    not have the correct statistics across processors. You can provide your own
305    parallel generator using PetscRandomRegister();
306 
307    If you create a PetscRandom() using PETSC_COMM_SELF on several processors then
308    the SAME random numbers will be generated on all those processors. Use PETSC_COMM_WORLD
309    or the appropriate parallel communicator to eliminate this issue.
310 
311    Use VecSetRandom() to set the elements of a vector to random numbers.
312 
313    Example of Usage:
314 .vb
315       PetscRandomCreate(PETSC_COMM_SELF,&r);
316       PetscRandomSetType(r,PETSCRAND48);
317       PetscRandomGetValue(r,&value1);
318       PetscRandomGetValueReal(r,&value2);
319       PetscRandomDestroy(&r);
320 .ve
321 
322    Concepts: random numbers^creating
323 
324 .seealso: PetscRandomSetType(), PetscRandomGetValue(), PetscRandomGetValueReal(), PetscRandomSetInterval(),
325           PetscRandomDestroy(), VecSetRandom(), PetscRandomType
326 @*/
327 
328 PetscErrorCode  PetscRandomCreate(MPI_Comm comm,PetscRandom *r)
329 {
330   PetscRandom    rr;
331   PetscErrorCode ierr;
332   PetscMPIInt    rank;
333 
334   PetscFunctionBegin;
335   PetscValidPointer(r,3);
336   *r = NULL;
337   ierr = PetscRandomInitializePackage();CHKERRQ(ierr);
338 
339   ierr = PetscHeaderCreate(rr,_p_PetscRandom,struct _PetscRandomOps,PETSC_RANDOM_CLASSID,"PetscRandom","Random number generator","Sys",comm,PetscRandomDestroy,0);CHKERRQ(ierr);
340 
341   ierr = MPI_Comm_rank(comm,&rank);CHKERRQ(ierr);
342 
343   rr->data  = NULL;
344   rr->low   = 0.0;
345   rr->width = 1.0;
346   rr->iset  = PETSC_FALSE;
347   rr->seed  = 0x12345678 + 76543*rank;
348   *r = rr;
349   PetscFunctionReturn(0);
350 }
351 
352 #undef __FUNCT__
353 #define __FUNCT__ "PetscRandomSeed"
354 /*@
355    PetscRandomSeed - Seed the generator.
356 
357    Not collective
358 
359    Input Parameters:
360 .  r - The random number generator context
361 
362    Level: intermediate
363 
364    Usage:
365       PetscRandomSetSeed(r,a positive integer);
366       PetscRandomSeed(r);  PetscRandomGetValue() will now start with the new seed.
367 
368       PetscRandomSeed(r) without a call to PetscRandomSetSeed() re-initializes
369         the seed. The random numbers generated will be the same as before.
370 
371    Concepts: random numbers^seed
372 
373 .seealso: PetscRandomCreate(), PetscRandomGetSeed(), PetscRandomSetSeed()
374 @*/
375 PetscErrorCode  PetscRandomSeed(PetscRandom r)
376 {
377   PetscErrorCode ierr;
378 
379   PetscFunctionBegin;
380   PetscValidHeaderSpecific(r,PETSC_RANDOM_CLASSID,1);
381   PetscValidType(r,1);
382 
383   ierr = (*r->ops->seed)(r);CHKERRQ(ierr);
384   ierr = PetscObjectStateIncrease((PetscObject)r);CHKERRQ(ierr);
385   PetscFunctionReturn(0);
386 }
387 
388