xref: /petsc/src/sys/classes/random/interface/randomc.c (revision b3bd51de3a0e6305f514f78fecfad5e8b5d71828)
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     defaultType = PETSCRANDER48;
142   }
143 
144   ierr = PetscRandomRegisterAll();CHKERRQ(ierr);
145   ierr = PetscOptionsFList("-random_type","PetscRandom type","PetscRandomSetType",PetscRandomList,defaultType,typeName,256,&opt);CHKERRQ(ierr);
146   if (opt) {
147     ierr = PetscRandomSetType(rnd, typeName);CHKERRQ(ierr);
148   } else {
149     ierr = PetscRandomSetType(rnd, defaultType);CHKERRQ(ierr);
150   }
151   PetscFunctionReturn(0);
152 }
153 
154 #undef __FUNCT__
155 #define __FUNCT__ "PetscRandomSetFromOptions"
156 /*@
157   PetscRandomSetFromOptions - Configures the random number generator from the options database.
158 
159   Collective on PetscRandom
160 
161   Input Parameter:
162 . rnd - The random number generator context
163 
164   Options Database:
165 .  -random_seed <integer> - provide a seed to the random number generater
166 
167   Notes:  To see all options, run your program with the -help option.
168           Must be called after PetscRandomCreate() but before the rnd is used.
169 
170   Level: beginner
171 
172 .keywords: PetscRandom, set, options, database
173 .seealso: PetscRandomCreate(), PetscRandomSetType()
174 @*/
175 PetscErrorCode  PetscRandomSetFromOptions(PetscRandom rnd)
176 {
177   PetscErrorCode ierr;
178   PetscBool      set;
179   PetscInt       seed;
180 
181   PetscFunctionBegin;
182   PetscValidHeaderSpecific(rnd,PETSC_RANDOM_CLASSID,1);
183 
184   ierr = PetscObjectOptionsBegin((PetscObject)rnd);CHKERRQ(ierr);
185 
186   /* Handle PetscRandom type options */
187   ierr = PetscRandomSetTypeFromOptions_Private(PetscOptionsObject,rnd);CHKERRQ(ierr);
188 
189   /* Handle specific random generator's options */
190   if (rnd->ops->setfromoptions) {
191     ierr = (*rnd->ops->setfromoptions)(PetscOptionsObject,rnd);CHKERRQ(ierr);
192   }
193   ierr = PetscOptionsInt("-random_seed","Seed to use to generate random numbers","PetscRandomSetSeed",0,&seed,&set);CHKERRQ(ierr);
194   if (set) {
195     ierr = PetscRandomSetSeed(rnd,(unsigned long int)seed);CHKERRQ(ierr);
196     ierr = PetscRandomSeed(rnd);CHKERRQ(ierr);
197   }
198   ierr = PetscOptionsEnd();CHKERRQ(ierr);
199   ierr = PetscRandomViewFromOptions(rnd,NULL, "-random_view");CHKERRQ(ierr);
200   PetscFunctionReturn(0);
201 }
202 
203 #if defined(PETSC_HAVE_SAWS)
204 #include <petscviewersaws.h>
205 #endif
206 #undef __FUNCT__
207 #define __FUNCT__ "PetscRandomView"
208 /*@C
209    PetscRandomView - Views a random number generator object.
210 
211    Collective on PetscRandom
212 
213    Input Parameters:
214 +  rnd - The random number generator context
215 -  viewer - an optional visualization context
216 
217    Notes:
218    The available visualization contexts include
219 +     PETSC_VIEWER_STDOUT_SELF - standard output (default)
220 -     PETSC_VIEWER_STDOUT_WORLD - synchronized standard
221          output where only the first processor opens
222          the file.  All other processors send their
223          data to the first processor to print.
224 
225    You can change the format the vector is printed using the
226    option PetscViewerSetFormat().
227 
228    Level: beginner
229 
230 .seealso:  PetscRealView(), PetscScalarView(), PetscIntView()
231 @*/
232 PetscErrorCode  PetscRandomView(PetscRandom rnd,PetscViewer viewer)
233 {
234   PetscErrorCode ierr;
235   PetscBool      iascii;
236 #if defined(PETSC_HAVE_SAWS)
237   PetscBool      issaws;
238 #endif
239 
240   PetscFunctionBegin;
241   PetscValidHeaderSpecific(rnd,PETSC_RANDOM_CLASSID,1);
242   PetscValidType(rnd,1);
243   if (!viewer) {
244     ierr = PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)rnd),&viewer);CHKERRQ(ierr);
245   }
246   PetscValidHeaderSpecific(viewer,PETSC_VIEWER_CLASSID,2);
247   PetscCheckSameComm(rnd,1,viewer,2);
248   ierr = PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);CHKERRQ(ierr);
249 #if defined(PETSC_HAVE_SAWS)
250   ierr = PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERSAWS,&issaws);CHKERRQ(ierr);
251 #endif
252   if (iascii) {
253     PetscMPIInt rank;
254     ierr = PetscObjectPrintClassNamePrefixType((PetscObject)rnd,viewer);CHKERRQ(ierr);
255     ierr = MPI_Comm_rank(PetscObjectComm((PetscObject)rnd),&rank);CHKERRQ(ierr);
256     ierr = PetscViewerASCIIPushSynchronized(viewer);CHKERRQ(ierr);
257     ierr = PetscViewerASCIISynchronizedPrintf(viewer,"[%d] Random type %s, seed %D\n",rank,((PetscObject)rnd)->type_name,rnd->seed);CHKERRQ(ierr);
258     ierr = PetscViewerFlush(viewer);CHKERRQ(ierr);
259     ierr = PetscViewerASCIIPopSynchronized(viewer);CHKERRQ(ierr);
260 #if defined(PETSC_HAVE_SAWS)
261   } else if (issaws) {
262     PetscMPIInt rank;
263     const char  *name;
264 
265     ierr = PetscObjectGetName((PetscObject)rnd,&name);CHKERRQ(ierr);
266     ierr = MPI_Comm_rank(PETSC_COMM_WORLD,&rank);CHKERRQ(ierr);
267     if (!((PetscObject)rnd)->amsmem && !rank) {
268       char       dir[1024];
269 
270       ierr = PetscObjectViewSAWs((PetscObject)rnd,viewer);CHKERRQ(ierr);
271       ierr = PetscSNPrintf(dir,1024,"/PETSc/Objects/%s/Low",name);CHKERRQ(ierr);
272       PetscStackCallSAWs(SAWs_Register,(dir,&rnd->low,1,SAWs_READ,SAWs_DOUBLE));
273     }
274 #endif
275   }
276   PetscFunctionReturn(0);
277 }
278 
279 #undef __FUNCT__
280 #define __FUNCT__ "PetscRandomCreate"
281 /*@
282    PetscRandomCreate - Creates a context for generating random numbers,
283    and initializes the random-number generator.
284 
285    Collective on MPI_Comm
286 
287    Input Parameters:
288 +  comm - MPI communicator
289 
290    Output Parameter:
291 .  r  - the random number generator context
292 
293    Level: intermediate
294 
295    Notes:
296    The random type has to be set by PetscRandomSetType().
297 
298    This is only a primative "parallel" random number generator, it should NOT
299    be used for sophisticated parallel Monte Carlo methods since it will very likely
300    not have the correct statistics across processors. You can provide your own
301    parallel generator using PetscRandomRegister();
302 
303    If you create a PetscRandom() using PETSC_COMM_SELF on several processors then
304    the SAME random numbers will be generated on all those processors. Use PETSC_COMM_WORLD
305    or the appropriate parallel communicator to eliminate this issue.
306 
307    Use VecSetRandom() to set the elements of a vector to random numbers.
308 
309    Example of Usage:
310 .vb
311       PetscRandomCreate(PETSC_COMM_SELF,&r);
312       PetscRandomSetType(r,PETSCRAND48);
313       PetscRandomGetValue(r,&value1);
314       PetscRandomGetValueReal(r,&value2);
315       PetscRandomDestroy(&r);
316 .ve
317 
318    Concepts: random numbers^creating
319 
320 .seealso: PetscRandomSetType(), PetscRandomGetValue(), PetscRandomGetValueReal(), PetscRandomSetInterval(),
321           PetscRandomDestroy(), VecSetRandom(), PetscRandomType
322 @*/
323 
324 PetscErrorCode  PetscRandomCreate(MPI_Comm comm,PetscRandom *r)
325 {
326   PetscRandom    rr;
327   PetscErrorCode ierr;
328   PetscMPIInt    rank;
329 
330   PetscFunctionBegin;
331   PetscValidPointer(r,3);
332   *r = NULL;
333   ierr = PetscRandomInitializePackage();CHKERRQ(ierr);
334 
335   ierr = PetscHeaderCreate(rr,PETSC_RANDOM_CLASSID,"PetscRandom","Random number generator","Sys",comm,PetscRandomDestroy,NULL);CHKERRQ(ierr);
336 
337   ierr = MPI_Comm_rank(comm,&rank);CHKERRQ(ierr);
338 
339   rr->data  = NULL;
340   rr->low   = 0.0;
341   rr->width = 1.0;
342   rr->iset  = PETSC_FALSE;
343   rr->seed  = 0x12345678 + 76543*rank;
344   *r = rr;
345   PetscFunctionReturn(0);
346 }
347 
348 #undef __FUNCT__
349 #define __FUNCT__ "PetscRandomSeed"
350 /*@
351    PetscRandomSeed - Seed the generator.
352 
353    Not collective
354 
355    Input Parameters:
356 .  r - The random number generator context
357 
358    Level: intermediate
359 
360    Usage:
361       PetscRandomSetSeed(r,a positive integer);
362       PetscRandomSeed(r);  PetscRandomGetValue() will now start with the new seed.
363 
364       PetscRandomSeed(r) without a call to PetscRandomSetSeed() re-initializes
365         the seed. The random numbers generated will be the same as before.
366 
367    Concepts: random numbers^seed
368 
369 .seealso: PetscRandomCreate(), PetscRandomGetSeed(), PetscRandomSetSeed()
370 @*/
371 PetscErrorCode  PetscRandomSeed(PetscRandom r)
372 {
373   PetscErrorCode ierr;
374 
375   PetscFunctionBegin;
376   PetscValidHeaderSpecific(r,PETSC_RANDOM_CLASSID,1);
377   PetscValidType(r,1);
378 
379   ierr = (*r->ops->seed)(r);CHKERRQ(ierr);
380   ierr = PetscObjectStateIncrease((PetscObject)r);CHKERRQ(ierr);
381   PetscFunctionReturn(0);
382 }
383 
384