Actual source code: randomc.c
petsc-3.10.5 2019-03-28
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:
7: Numbers can be shuffled for use as random tuples
9: Multiple random number generators may be used
11: We are still not sure what interface we want here. There should be
12: one to reinitialize and set the seed.
13: */
15: #include <../src/sys/classes/random/randomimpl.h>
16: #include <petscviewer.h>
18: /* Logging support */
19: PetscClassId PETSC_RANDOM_CLASSID;
21: /*@
22: PetscRandomDestroy - Destroys a context that has been formed by
23: PetscRandomCreate().
25: Collective on PetscRandom
27: Intput Parameter:
28: . r - the random number generator context
30: Level: intermediate
32: .seealso: PetscRandomGetValue(), PetscRandomCreate(), VecSetRandom()
33: @*/
34: PetscErrorCode PetscRandomDestroy(PetscRandom *r)
35: {
39: if (!*r) return(0);
41: if (--((PetscObject)(*r))->refct > 0) {*r = 0; return(0);}
42: if ((*r)->ops->destroy) {
43: (*(*r)->ops->destroy)(*r);
44: }
45: PetscHeaderDestroy(r);
46: return(0);
47: }
50: /*@C
51: PetscRandomGetSeed - Gets the random seed.
53: Not collective
55: Input Parameters:
56: . r - The random number generator context
58: Output Parameter:
59: . seed - The random seed
61: Level: intermediate
63: Concepts: random numbers^seed
65: .seealso: PetscRandomCreate(), PetscRandomSetSeed(), PetscRandomSeed()
66: @*/
67: PetscErrorCode PetscRandomGetSeed(PetscRandom r,unsigned long *seed)
68: {
71: if (seed) {
73: *seed = r->seed;
74: }
75: return(0);
76: }
78: /*@C
79: PetscRandomSetSeed - Sets the random seed. You MUST call PetscRandomSeed() after this call to have the new seed used.
81: Not collective
83: Input Parameters:
84: + r - The random number generator context
85: - seed - The random seed
87: Level: intermediate
89: Usage:
90: PetscRandomSetSeed(r,a positive integer);
91: PetscRandomSeed(r); PetscRandomGetValue() will now start with the new seed.
93: PetscRandomSeed(r) without a call to PetscRandomSetSeed() re-initializes
94: the seed. The random numbers generated will be the same as before.
96: Concepts: random numbers^seed
98: .seealso: PetscRandomCreate(), PetscRandomGetSeed(), PetscRandomSeed()
99: @*/
100: PetscErrorCode PetscRandomSetSeed(PetscRandom r,unsigned long seed)
101: {
106: r->seed = seed;
107: PetscInfo1(NULL,"Setting seed to %d\n",(int)seed);
108: return(0);
109: }
111: /* ------------------------------------------------------------------- */
112: /*
113: PetscRandomSetTypeFromOptions_Private - Sets the type of random generator from user options. Defaults to type PETSCRAND48 or PETSCRAND.
115: Collective on PetscRandom
117: Input Parameter:
118: . rnd - The random number generator context
120: Level: intermediate
122: .keywords: PetscRandom, set, options, database, type
123: .seealso: PetscRandomSetFromOptions(), PetscRandomSetType()
124: */
125: static PetscErrorCode PetscRandomSetTypeFromOptions_Private(PetscOptionItems *PetscOptionsObject,PetscRandom rnd)
126: {
127: PetscBool opt;
128: const char *defaultType;
129: char typeName[256];
133: if (((PetscObject)rnd)->type_name) {
134: defaultType = ((PetscObject)rnd)->type_name;
135: } else {
136: defaultType = PETSCRANDER48;
137: }
139: PetscRandomRegisterAll();
140: PetscOptionsFList("-random_type","PetscRandom type","PetscRandomSetType",PetscRandomList,defaultType,typeName,256,&opt);
141: if (opt) {
142: PetscRandomSetType(rnd, typeName);
143: } else {
144: PetscRandomSetType(rnd, defaultType);
145: }
146: return(0);
147: }
149: /*@
150: PetscRandomSetFromOptions - Configures the random number generator from the options database.
152: Collective on PetscRandom
154: Input Parameter:
155: . rnd - The random number generator context
157: Options Database:
158: + -random_seed <integer> - provide a seed to the random number generater
159: - -random_no_imaginary_part - makes the imaginary part of the random number zero, this is useful when you want the
160: same code to produce the same result when run with real numbers or complex numbers for regression testing purposes
162: Notes:
163: To see all options, run your program with the -help option.
164: Must be called after PetscRandomCreate() but before the rnd is used.
166: Level: beginner
168: .keywords: PetscRandom, set, options, database
169: .seealso: PetscRandomCreate(), PetscRandomSetType()
170: @*/
171: PetscErrorCode PetscRandomSetFromOptions(PetscRandom rnd)
172: {
174: PetscBool set,noimaginary = PETSC_FALSE;
175: PetscInt seed;
180: PetscObjectOptionsBegin((PetscObject)rnd);
182: /* Handle PetscRandom type options */
183: PetscRandomSetTypeFromOptions_Private(PetscOptionsObject,rnd);
185: /* Handle specific random generator's options */
186: if (rnd->ops->setfromoptions) {
187: (*rnd->ops->setfromoptions)(PetscOptionsObject,rnd);
188: }
189: PetscOptionsInt("-random_seed","Seed to use to generate random numbers","PetscRandomSetSeed",0,&seed,&set);
190: if (set) {
191: PetscRandomSetSeed(rnd,(unsigned long int)seed);
192: PetscRandomSeed(rnd);
193: }
194: PetscOptionsBool("-random_no_imaginary_part","The imaginary part of the random number will be zero","PetscRandomSetInterval",noimaginary,&noimaginary,&set);
195: #if defined(PETSC_HAVE_COMPLEX)
196: if (set) {
197: if (noimaginary) {
198: PetscScalar low,high;
199: PetscRandomGetInterval(rnd,&low,&high);
200: low = low - PetscImaginaryPart(low);
201: high = high - PetscImaginaryPart(high);
202: PetscRandomSetInterval(rnd,low,high);
203: }
204: }
205: #endif
206: PetscOptionsEnd();
207: PetscRandomViewFromOptions(rnd,NULL, "-random_view");
208: return(0);
209: }
211: #if defined(PETSC_HAVE_SAWS)
212: #include <petscviewersaws.h>
213: #endif
214: /*@C
215: PetscRandomView - Views a random number generator object.
217: Collective on PetscRandom
219: Input Parameters:
220: + rnd - The random number generator context
221: - viewer - an optional visualization context
223: Notes:
224: The available visualization contexts include
225: + PETSC_VIEWER_STDOUT_SELF - standard output (default)
226: - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
227: output where only the first processor opens
228: the file. All other processors send their
229: data to the first processor to print.
231: You can change the format the vector is printed using the
232: option PetscViewerPushFormat().
234: Level: beginner
236: .seealso: PetscRealView(), PetscScalarView(), PetscIntView()
237: @*/
238: PetscErrorCode PetscRandomView(PetscRandom rnd,PetscViewer viewer)
239: {
241: PetscBool iascii;
242: #if defined(PETSC_HAVE_SAWS)
243: PetscBool issaws;
244: #endif
249: if (!viewer) {
250: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)rnd),&viewer);
251: }
254: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
255: #if defined(PETSC_HAVE_SAWS)
256: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERSAWS,&issaws);
257: #endif
258: if (iascii) {
259: PetscMPIInt rank;
260: PetscObjectPrintClassNamePrefixType((PetscObject)rnd,viewer);
261: MPI_Comm_rank(PetscObjectComm((PetscObject)rnd),&rank);
262: PetscViewerASCIIPushSynchronized(viewer);
263: PetscViewerASCIISynchronizedPrintf(viewer,"[%d] Random type %s, seed %lu\n",rank,((PetscObject)rnd)->type_name,rnd->seed);
264: PetscViewerFlush(viewer);
265: PetscViewerASCIIPopSynchronized(viewer);
266: #if defined(PETSC_HAVE_SAWS)
267: } else if (issaws) {
268: PetscMPIInt rank;
269: const char *name;
271: PetscObjectGetName((PetscObject)rnd,&name);
272: MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
273: if (!((PetscObject)rnd)->amsmem && !rank) {
274: char dir[1024];
276: PetscObjectViewSAWs((PetscObject)rnd,viewer);
277: PetscSNPrintf(dir,1024,"/PETSc/Objects/%s/Low",name);
278: PetscStackCallSAWs(SAWs_Register,(dir,&rnd->low,1,SAWs_READ,SAWs_DOUBLE));
279: }
280: #endif
281: }
282: return(0);
283: }
285: /*@
286: PetscRandomCreate - Creates a context for generating random numbers,
287: and initializes the random-number generator.
289: Collective on MPI_Comm
291: Input Parameters:
292: . comm - MPI communicator
294: Output Parameter:
295: . r - the random number generator context
297: Level: intermediate
299: Notes:
300: The random type has to be set by PetscRandomSetType().
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();
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.
311: Use VecSetRandom() to set the elements of a vector to random numbers.
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
322: Concepts: random numbers^creating
324: .seealso: PetscRandomSetType(), PetscRandomGetValue(), PetscRandomGetValueReal(), PetscRandomSetInterval(),
325: PetscRandomDestroy(), VecSetRandom(), PetscRandomType
326: @*/
328: PetscErrorCode PetscRandomCreate(MPI_Comm comm,PetscRandom *r)
329: {
330: PetscRandom rr;
332: PetscMPIInt rank;
336: *r = NULL;
337: PetscRandomInitializePackage();
339: PetscHeaderCreate(rr,PETSC_RANDOM_CLASSID,"PetscRandom","Random number generator","Sys",comm,PetscRandomDestroy,PetscRandomView);
341: MPI_Comm_rank(comm,&rank);
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: PetscRandomSetType(rr,PETSCRANDER48);
349: *r = rr;
350: return(0);
351: }
353: /*@
354: PetscRandomSeed - Seed the generator.
356: Not collective
358: Input Parameters:
359: . r - The random number generator context
361: Level: intermediate
363: Usage:
364: PetscRandomSetSeed(r,a positive integer);
365: PetscRandomSeed(r); PetscRandomGetValue() will now start with the new seed.
367: PetscRandomSeed(r) without a call to PetscRandomSetSeed() re-initializes
368: the seed. The random numbers generated will be the same as before.
370: Concepts: random numbers^seed
372: .seealso: PetscRandomCreate(), PetscRandomGetSeed(), PetscRandomSetSeed()
373: @*/
374: PetscErrorCode PetscRandomSeed(PetscRandom r)
375: {
382: (*r->ops->seed)(r);
383: PetscObjectStateIncrease((PetscObject)r);
384: return(0);
385: }