Actual source code: randomc.c
petsc-3.7.5 2017-01-01
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> /*I "petscsys.h" I*/
16: #include <petscviewer.h>
18: /* Logging support */
19: PetscClassId PETSC_RANDOM_CLASSID;
23: /*@
24: PetscRandomDestroy - Destroys a context that has been formed by
25: PetscRandomCreate().
27: Collective on PetscRandom
29: Intput Parameter:
30: . r - the random number generator context
32: Level: intermediate
34: .seealso: PetscRandomGetValue(), PetscRandomCreate(), VecSetRandom()
35: @*/
36: PetscErrorCode PetscRandomDestroy(PetscRandom *r)
37: {
41: if (!*r) return(0);
43: if (--((PetscObject)(*r))->refct > 0) {*r = 0; return(0);}
44: if ((*r)->ops->destroy) {
45: (*(*r)->ops->destroy)(*r);
46: }
47: PetscHeaderDestroy(r);
48: return(0);
49: }
54: /*@
55: PetscRandomGetSeed - Gets the random seed.
57: Not collective
59: Input Parameters:
60: . r - The random number generator context
62: Output Parameter:
63: . seed - The random seed
65: Level: intermediate
67: Concepts: random numbers^seed
69: .seealso: PetscRandomCreate(), PetscRandomSetSeed(), PetscRandomSeed()
70: @*/
71: PetscErrorCode PetscRandomGetSeed(PetscRandom r,unsigned long *seed)
72: {
75: if (seed) {
77: *seed = r->seed;
78: }
79: return(0);
80: }
84: /*@
85: PetscRandomSetSeed - Sets the random seed. You MUST call PetscRandomSeed() after this call to have the new seed used.
87: Not collective
89: Input Parameters:
90: + r - The random number generator context
91: - seed - The random seed
93: Level: intermediate
95: Usage:
96: PetscRandomSetSeed(r,a positive integer);
97: PetscRandomSeed(r); PetscRandomGetValue() will now start with the new seed.
99: PetscRandomSeed(r) without a call to PetscRandomSetSeed() re-initializes
100: the seed. The random numbers generated will be the same as before.
102: Concepts: random numbers^seed
104: .seealso: PetscRandomCreate(), PetscRandomGetSeed(), PetscRandomSeed()
105: @*/
106: PetscErrorCode PetscRandomSetSeed(PetscRandom r,unsigned long seed)
107: {
112: r->seed = seed;
113: PetscInfo1(NULL,"Setting seed to %d\n",(int)seed);
114: return(0);
115: }
117: /* ------------------------------------------------------------------- */
120: /*
121: PetscRandomSetTypeFromOptions_Private - Sets the type of random generator from user options. Defaults to type PETSCRAND48 or PETSCRAND.
123: Collective on PetscRandom
125: Input Parameter:
126: . rnd - The random number generator context
128: Level: intermediate
130: .keywords: PetscRandom, set, options, database, type
131: .seealso: PetscRandomSetFromOptions(), PetscRandomSetType()
132: */
133: static PetscErrorCode PetscRandomSetTypeFromOptions_Private(PetscOptionItems *PetscOptionsObject,PetscRandom rnd)
134: {
135: PetscBool opt;
136: const char *defaultType;
137: char typeName[256];
141: if (((PetscObject)rnd)->type_name) {
142: defaultType = ((PetscObject)rnd)->type_name;
143: } else {
144: defaultType = PETSCRANDER48;
145: }
147: PetscRandomRegisterAll();
148: PetscOptionsFList("-random_type","PetscRandom type","PetscRandomSetType",PetscRandomList,defaultType,typeName,256,&opt);
149: if (opt) {
150: PetscRandomSetType(rnd, typeName);
151: } else {
152: PetscRandomSetType(rnd, defaultType);
153: }
154: return(0);
155: }
159: /*@
160: PetscRandomSetFromOptions - Configures the random number generator from the options database.
162: Collective on PetscRandom
164: Input Parameter:
165: . rnd - The random number generator context
167: Options Database:
168: + -random_seed <integer> - provide a seed to the random number generater
169: - -random_no_imaginary_part - makes the imaginary part of the random number zero, this is useful when you want the
170: same code to produce the same result when run with real numbers or complex numbers for regression testing purposes
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.
175: Level: beginner
177: .keywords: PetscRandom, set, options, database
178: .seealso: PetscRandomCreate(), PetscRandomSetType()
179: @*/
180: PetscErrorCode PetscRandomSetFromOptions(PetscRandom rnd)
181: {
183: PetscBool set,noimaginary = PETSC_FALSE;
184: PetscInt seed;
189: PetscObjectOptionsBegin((PetscObject)rnd);
191: /* Handle PetscRandom type options */
192: PetscRandomSetTypeFromOptions_Private(PetscOptionsObject,rnd);
194: /* Handle specific random generator's options */
195: if (rnd->ops->setfromoptions) {
196: (*rnd->ops->setfromoptions)(PetscOptionsObject,rnd);
197: }
198: PetscOptionsInt("-random_seed","Seed to use to generate random numbers","PetscRandomSetSeed",0,&seed,&set);
199: if (set) {
200: PetscRandomSetSeed(rnd,(unsigned long int)seed);
201: PetscRandomSeed(rnd);
202: }
203: PetscOptionsBool("-random_no_imaginary_part","The imaginary part of the random number will be zero","PetscRandomSetInterval",noimaginary,&noimaginary,&set);
204: #if defined(PETSC_HAVE_COMPLEX)
205: if (set) {
206: if (noimaginary) {
207: PetscScalar low,high;
208: PetscRandomGetInterval(rnd,&low,&high);
209: low = low - PetscImaginaryPart(low);
210: high = high - PetscImaginaryPart(high);
211: PetscRandomSetInterval(rnd,low,high);
212: }
213: }
214: #endif
215: PetscOptionsEnd();
216: PetscRandomViewFromOptions(rnd,NULL, "-random_view");
217: return(0);
218: }
220: #if defined(PETSC_HAVE_SAWS)
221: #include <petscviewersaws.h>
222: #endif
225: /*@C
226: PetscRandomView - Views a random number generator object.
228: Collective on PetscRandom
230: Input Parameters:
231: + rnd - The random number generator context
232: - viewer - an optional visualization context
234: Notes:
235: The available visualization contexts include
236: + PETSC_VIEWER_STDOUT_SELF - standard output (default)
237: - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
238: output where only the first processor opens
239: the file. All other processors send their
240: data to the first processor to print.
242: You can change the format the vector is printed using the
243: option PetscViewerPushFormat().
245: Level: beginner
247: .seealso: PetscRealView(), PetscScalarView(), PetscIntView()
248: @*/
249: PetscErrorCode PetscRandomView(PetscRandom rnd,PetscViewer viewer)
250: {
252: PetscBool iascii;
253: #if defined(PETSC_HAVE_SAWS)
254: PetscBool issaws;
255: #endif
260: if (!viewer) {
261: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)rnd),&viewer);
262: }
265: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
266: #if defined(PETSC_HAVE_SAWS)
267: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERSAWS,&issaws);
268: #endif
269: if (iascii) {
270: PetscMPIInt rank;
271: PetscObjectPrintClassNamePrefixType((PetscObject)rnd,viewer);
272: MPI_Comm_rank(PetscObjectComm((PetscObject)rnd),&rank);
273: PetscViewerASCIIPushSynchronized(viewer);
274: PetscViewerASCIISynchronizedPrintf(viewer,"[%d] Random type %s, seed %D\n",rank,((PetscObject)rnd)->type_name,rnd->seed);
275: PetscViewerFlush(viewer);
276: PetscViewerASCIIPopSynchronized(viewer);
277: #if defined(PETSC_HAVE_SAWS)
278: } else if (issaws) {
279: PetscMPIInt rank;
280: const char *name;
282: PetscObjectGetName((PetscObject)rnd,&name);
283: MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
284: if (!((PetscObject)rnd)->amsmem && !rank) {
285: char dir[1024];
287: PetscObjectViewSAWs((PetscObject)rnd,viewer);
288: PetscSNPrintf(dir,1024,"/PETSc/Objects/%s/Low",name);
289: PetscStackCallSAWs(SAWs_Register,(dir,&rnd->low,1,SAWs_READ,SAWs_DOUBLE));
290: }
291: #endif
292: }
293: return(0);
294: }
298: /*@
299: PetscRandomCreate - Creates a context for generating random numbers,
300: and initializes the random-number generator.
302: Collective on MPI_Comm
304: Input Parameters:
305: + comm - MPI communicator
307: Output Parameter:
308: . r - the random number generator context
310: Level: intermediate
312: Notes:
313: The random type has to be set by PetscRandomSetType().
315: This is only a primative "parallel" random number generator, it should NOT
316: be used for sophisticated parallel Monte Carlo methods since it will very likely
317: not have the correct statistics across processors. You can provide your own
318: parallel generator using PetscRandomRegister();
320: If you create a PetscRandom() using PETSC_COMM_SELF on several processors then
321: the SAME random numbers will be generated on all those processors. Use PETSC_COMM_WORLD
322: or the appropriate parallel communicator to eliminate this issue.
324: Use VecSetRandom() to set the elements of a vector to random numbers.
326: Example of Usage:
327: .vb
328: PetscRandomCreate(PETSC_COMM_SELF,&r);
329: PetscRandomSetType(r,PETSCRAND48);
330: PetscRandomGetValue(r,&value1);
331: PetscRandomGetValueReal(r,&value2);
332: PetscRandomDestroy(&r);
333: .ve
335: Concepts: random numbers^creating
337: .seealso: PetscRandomSetType(), PetscRandomGetValue(), PetscRandomGetValueReal(), PetscRandomSetInterval(),
338: PetscRandomDestroy(), VecSetRandom(), PetscRandomType
339: @*/
341: PetscErrorCode PetscRandomCreate(MPI_Comm comm,PetscRandom *r)
342: {
343: PetscRandom rr;
345: PetscMPIInt rank;
349: *r = NULL;
350: PetscRandomInitializePackage();
352: PetscHeaderCreate(rr,PETSC_RANDOM_CLASSID,"PetscRandom","Random number generator","Sys",comm,PetscRandomDestroy,PetscRandomView);
354: MPI_Comm_rank(comm,&rank);
356: rr->data = NULL;
357: rr->low = 0.0;
358: rr->width = 1.0;
359: rr->iset = PETSC_FALSE;
360: rr->seed = 0x12345678 + 76543*rank;
361: PetscRandomSetType(rr,PETSCRANDER48);
362: *r = rr;
363: return(0);
364: }
368: /*@
369: PetscRandomSeed - Seed the generator.
371: Not collective
373: Input Parameters:
374: . r - The random number generator context
376: Level: intermediate
378: Usage:
379: PetscRandomSetSeed(r,a positive integer);
380: PetscRandomSeed(r); PetscRandomGetValue() will now start with the new seed.
382: PetscRandomSeed(r) without a call to PetscRandomSetSeed() re-initializes
383: the seed. The random numbers generated will be the same as before.
385: Concepts: random numbers^seed
387: .seealso: PetscRandomCreate(), PetscRandomGetSeed(), PetscRandomSetSeed()
388: @*/
389: PetscErrorCode PetscRandomSeed(PetscRandom r)
390: {
397: (*r->ops->seed)(r);
398: PetscObjectStateIncrease((PetscObject)r);
399: return(0);
400: }