Actual source code: random.c
1: #define PETSC_DLL
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/random/randomimpl.h
16: #if defined (PETSC_HAVE_STDLIB_H)
17: #include <stdlib.h>
18: #endif
20: /* Logging support */
21: PetscCookie PETSC_RANDOM_COOKIE;
25: /*@
26: PetscRandomDestroy - Destroys a context that has been formed by
27: PetscRandomCreate().
29: Collective on PetscRandom
31: Intput Parameter:
32: . r - the random number generator context
34: Level: intermediate
36: .seealso: PetscRandomGetValue(), PetscRandomCreate(), VecSetRandom()
37: @*/
38: PetscErrorCode PetscRandomDestroy(PetscRandom r)
39: {
43: if (--((PetscObject)r)->refct > 0) return(0);
44: PetscHeaderDestroy(r);
45: return(0);
46: }
50: /*@C
51: PetscRandomGetInterval - Gets the interval over which the random numbers
52: will be randomly distributed. By default, this interval is [0,1).
54: Not collective
56: Input Parameters:
57: . r - the random number generator context
59: Output Parameters:
60: + low - The lower bound of the interval
61: - high - The upper bound of the interval
63: Level: intermediate
65: Concepts: random numbers^range
67: .seealso: PetscRandomCreate(), PetscRandomSetInterval()
68: @*/
69: PetscErrorCode PetscRandomGetInterval(PetscRandom r,PetscScalar *low,PetscScalar *high)
70: {
73: if (low) {
75: *low = r->low;
76: }
77: if (high) {
79: *high = r->low+r->width;
80: }
81: return(0);
82: }
86: /*@
87: PetscRandomSetInterval - Sets the interval over which the random numbers
88: will be randomly distributed. By default, this interval is [0,1).
90: Not collective
92: Input Parameters:
93: + r - the random number generator context
94: . low - The lower bound of the interval
95: - high - The upper bound of the interval
97: Level: intermediate
99: Notes: for complex numbers either the real part or the imaginary part of high must be greater than its low part; or both of them can be greater.
100: If the real or imaginary part of low and high are the same then that value is always returned in the real or imaginary part.
102: Concepts: random numbers^range
104: .seealso: PetscRandomCreate(), PetscRandomGetInterval()
105: @*/
106: PetscErrorCode PetscRandomSetInterval(PetscRandom r,PetscScalar low,PetscScalar high)
107: {
110: #if defined(PETSC_USE_COMPLEX)
111: if (PetscRealPart(low) > PetscRealPart(high)) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"only low < high");
112: if (PetscImaginaryPart(low) > PetscImaginaryPart(high)) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"only low < high");
113: #else
114: if (low >= high) SETERRQ2(PETSC_ERR_ARG_OUTOFRANGE,"only low < high: Instead %G %G",low,high);
115: #endif
116: r->low = low;
117: r->width = high-low;
118: r->iset = PETSC_TRUE;
119: return(0);
120: }
124: /*@C
125: PetscRandomGetSeed - Gets the random seed.
127: Not collective
129: Input Parameters:
130: . r - The random number generator context
132: Output Parameter:
133: . seed - The random seed
135: Level: intermediate
137: Concepts: random numbers^seed
139: .seealso: PetscRandomCreate(), PetscRandomSetSeed(), PetscRandomSeed()
140: @*/
141: PetscErrorCode PetscRandomGetSeed(PetscRandom r,unsigned long *seed)
142: {
145: if (seed) {
147: *seed = r->seed;
148: }
149: return(0);
150: }
154: /*@C
155: PetscRandomSetSeed - Sets the random seed.
157: Not collective
159: Input Parameters:
160: + r - The random number generator context
161: - seed - The random seed
163: Level: intermediate
165: Usage:
166: PetscRandomSetSeed(r,a positive integer);
167: PetscRandomSeed(r); PetscRandomGetValue() will now start with the new seed.
169: PetscRandomSeed(r) without a call to PetscRandomSetSeed() re-initializes
170: the seed. The random numbers generated will be the same as before.
172: Concepts: random numbers^seed
174: .seealso: PetscRandomCreate(), PetscRandomGetSeed(), PetscRandomSeed()
175: @*/
176: PetscErrorCode PetscRandomSetSeed(PetscRandom r,unsigned long seed)
177: {
180: r->seed = seed;
181: return(0);
182: }
184: /* ------------------------------------------------------------------- */
187: /*
188: PetscRandomSetTypeFromOptions_Private - Sets the type of random generator from user options. Defaults to type PETSCRAND48 or PETSCRAND.
190: Collective on PetscRandom
192: Input Parameter:
193: . rnd - The random number generator context
195: Level: intermediate
197: .keywords: PetscRandom, set, options, database, type
198: .seealso: PetscRandomSetFromOptions(), PetscRandomSetType()
199: */
200: static PetscErrorCode PetscRandomSetTypeFromOptions_Private(PetscRandom rnd)
201: {
202: PetscTruth opt;
203: const char *defaultType;
204: char typeName[256];
208: if (((PetscObject)rnd)->type_name) {
209: defaultType = ((PetscObject)rnd)->type_name;
210: } else {
211: #if defined(PETSC_HAVE_DRAND48)
212: defaultType = PETSCRAND48;
213: #elif defined(PETSC_HAVE_RAND)
214: defaultType = PETSCRAND;
215: #endif
216: }
218: if (!PetscRandomRegisterAllCalled) {PetscRandomRegisterAll(PETSC_NULL);}
219: PetscOptionsList("-random_type","PetscRandom type","PetscRandomSetType",PetscRandomList,defaultType,typeName,256,&opt);
220: if (opt) {
221: PetscRandomSetType(rnd, typeName);
222: } else {
223: PetscRandomSetType(rnd, defaultType);
224: }
225: return(0);
226: }
230: /*@
231: PetscRandomSetFromOptions - Configures the random number generator from the options database.
233: Collective on PetscRandom
235: Input Parameter:
236: . rnd - The random number generator context
238: Notes: To see all options, run your program with the -help option, or consult the users manual.
239: Must be called after PetscRandomCreate() but before the rnd is used.
241: Level: beginner
243: .keywords: PetscRandom, set, options, database
244: .seealso: PetscRandomCreate(), PetscRandomSetType()
245: @*/
246: PetscErrorCode PetscRandomSetFromOptions(PetscRandom rnd)
247: {
253: PetscOptionsBegin(((PetscObject)rnd)->comm, ((PetscObject)rnd)->prefix, "PetscRandom options", "PetscRandom");
255: /* Handle PetscRandom type options */
256: PetscRandomSetTypeFromOptions_Private(rnd);
258: /* Handle specific random generator's options */
259: if (rnd->ops->setfromoptions) {
260: (*rnd->ops->setfromoptions)(rnd);
261: }
262: PetscOptionsEnd();
263: PetscRandomViewFromOptions(rnd, ((PetscObject)rnd)->name);
264: return(0);
265: }
269: /*@C
270: PetscRandomView - Views a random number generator object.
272: Collective on PetscRandom
274: Input Parameters:
275: + rnd - The random number generator context
276: - viewer - an optional visualization context
278: Notes:
279: The available visualization contexts include
280: + PETSC_VIEWER_STDOUT_SELF - standard output (default)
281: - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
282: output where only the first processor opens
283: the file. All other processors send their
284: data to the first processor to print.
286: You can change the format the vector is printed using the
287: option PetscViewerSetFormat().
289: Level: beginner
291: .seealso: PetscRealView(), PetscScalarView(), PetscIntView()
292: @*/
293: PetscErrorCode PetscRandomView(PetscRandom rnd,PetscViewer viewer)
294: {
295: PetscErrorCode ierr;
296: PetscTruth iascii;
301: if (!viewer) {
302: PetscViewerASCIIGetStdout(((PetscObject)rnd)->comm,&viewer);
303: }
306: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&iascii);
307: if (iascii) {
308: PetscMPIInt rank;
309: MPI_Comm_rank(((PetscObject)rnd)->comm,&rank);
310: PetscViewerASCIISynchronizedPrintf(viewer,"[%D] Random type %s, seed %D\n",rank,((PetscObject)rnd)->type_name,rnd->seed);
311: PetscViewerFlush(viewer);
312: } else {
313: const char *tname;
314: PetscObjectGetName((PetscObject)viewer,&tname);
315: SETERRQ1(PETSC_ERR_SUP,"Viewer type %s not supported for this object",tname);
316: }
317: return(0);
318: }
320: #undef __FUNCT__
322: /*@
323: PetscRandomViewFromOptions - This function visualizes the type and the seed of the generated random numbers based upon user options.
325: Collective on PetscRandom
327: Input Parameters:
328: . rnd - The random number generator context
329: . title - The title
331: Level: intermediate
333: .keywords: PetscRandom, view, options, database
334: .seealso: PetscRandomSetFromOptions()
335: @*/
336: PetscErrorCode PetscRandomViewFromOptions(PetscRandom rnd, char *title)
337: {
338: PetscTruth opt = PETSC_FALSE;
339: PetscViewer viewer;
340: char typeName[1024];
341: char fileName[PETSC_MAX_PATH_LEN];
342: size_t len;
343:
347: PetscOptionsGetTruth(((PetscObject)rnd)->prefix, "-random_view", &opt,PETSC_NULL);
348: if (opt) {
349: PetscOptionsGetString(((PetscObject)rnd)->prefix, "-random_view", typeName, 1024, &opt);
350: PetscStrlen(typeName, &len);
351: if (len > 0) {
352: PetscViewerCreate(((PetscObject)rnd)->comm, &viewer);
353: PetscViewerSetType(viewer, typeName);
354: PetscOptionsGetString(((PetscObject)rnd)->prefix, "-random_view_file", fileName, 1024, &opt);
355: if (opt) {
356: PetscViewerFileSetName(viewer, fileName);
357: } else {
358: PetscViewerFileSetName(viewer, ((PetscObject)rnd)->name);
359: }
360: PetscRandomView(rnd, viewer);
361: PetscViewerFlush(viewer);
362: PetscViewerDestroy(viewer);
363: } else {
364: PetscViewer viewer;
365: PetscViewerASCIIGetStdout(((PetscObject)rnd)->comm,&viewer);
366: PetscRandomView(rnd, viewer);
367: }
368: }
369: return(0);
370: }
375: /*@
376: PetscRandomCreate - Creates a context for generating random numbers,
377: and initializes the random-number generator.
379: Collective on MPI_Comm
381: Input Parameters:
382: + comm - MPI communicator
384: Output Parameter:
385: . r - the random number generator context
387: Level: intermediate
389: Notes:
390: The random type has to be set by PetscRandomSetType().
392: This is only a primative "parallel" random number generator, it should NOT
393: be used for sophisticated parallel Monte Carlo methods since it will very likely
394: not have the correct statistics across processors. You can provide your own
395: parallel generator using PetscRandomRegister();
397: If you create a PetscRandom() using PETSC_COMM_SELF on several processors then
398: the SAME random numbers will be generated on all those processors. Use PETSC_COMM_WORLD
399: or the appropriate parallel communicator to eliminate this issue.
401: Use VecSetRandom() to set the elements of a vector to random numbers.
403: Example of Usage:
404: .vb
405: PetscRandomCreate(PETSC_COMM_SELF,&r);
406: PetscRandomSetType(r,PETSCRAND48);
407: PetscRandomGetValue(r,&value1);
408: PetscRandomGetValueReal(r,&value2);
409: PetscRandomDestroy(r);
410: .ve
412: Concepts: random numbers^creating
414: .seealso: PetscRandomSetType(), PetscRandomGetValue(), PetscRandomGetValueReal(), PetscRandomSetInterval(),
415: PetscRandomDestroy(), VecSetRandom(), PetscRandomType
416: @*/
418: PetscErrorCode PetscRandomCreate(MPI_Comm comm,PetscRandom *r)
419: {
420: PetscRandom rr;
422: PetscMPIInt rank;
426: *r = PETSC_NULL;
427: #ifndef PETSC_USE_DYNAMIC_LIBRARIES
428: PetscRandomInitializePackage(PETSC_NULL);
429: #endif
431: PetscHeaderCreate(rr,_p_PetscRandom,struct _PetscRandomOps,PETSC_RANDOM_COOKIE,-1,"PetscRandom",comm,PetscRandomDestroy,0);
433: MPI_Comm_rank(comm,&rank);
434: rr->data = PETSC_NULL;
435: rr->low = 0.0;
436: rr->width = 1.0;
437: rr->iset = PETSC_FALSE;
438: rr->seed = 0x12345678 + 76543*rank;
439: *r = rr;
440: return(0);
441: }
445: /*@C
446: PetscRandomSeed - Seed the generator.
448: Not collective
450: Input Parameters:
451: . r - The random number generator context
453: Level: intermediate
455: Usage:
456: PetscRandomSetSeed(r,a positive integer);
457: PetscRandomSeed(r); PetscRandomGetValue() will now start with the new seed.
459: PetscRandomSeed(r) without a call to PetscRandomSetSeed() re-initializes
460: the seed. The random numbers generated will be the same as before.
462: Concepts: random numbers^seed
464: .seealso: PetscRandomCreate(), PetscRandomGetSeed(), PetscRandomSetSeed()
465: @*/
466: PetscErrorCode PetscRandomSeed(PetscRandom r)
467: {
474: (*r->ops->seed)(r);
475: PetscObjectStateIncrease((PetscObject)r);
476: return(0);
477: }
481: /*@
482: PetscRandomGetValue - Generates a random number. Call this after first calling
483: PetscRandomCreate().
485: Not Collective
487: Intput Parameter:
488: . r - the random number generator context
490: Output Parameter:
491: . val - the value
493: Level: intermediate
495: Notes:
496: Use VecSetRandom() to set the elements of a vector to random numbers.
498: When PETSc is compiled for complex numbers this returns a complex number with random real and complex parts.
499: Use PetscGetValueReal() to get a random real number.
501: To get a complex number with only a random real part, first call PetscRandomSetInterval() with a equal
502: low and high imaginary part. Similarly to get a complex number with only a random imaginary part call
503: PetscRandomSetInterval() with a equal low and high real part.
505: Example of Usage:
506: .vb
507: PetscRandomCreate(PETSC_COMM_WORLD,&r);
508: PetscRandomGetValue(r,&value1);
509: PetscRandomGetValue(r,&value2);
510: PetscRandomGetValue(r,&value3);
511: PetscRandomDestroy(r);
512: .ve
514: Concepts: random numbers^getting
516: .seealso: PetscRandomCreate(), PetscRandomDestroy(), VecSetRandom(), PetscRandomGetValueReal()
517: @*/
518: PetscErrorCode PetscRandomGetValue(PetscRandom r,PetscScalar *val)
519: {
527: (*r->ops->getvalue)(r,val);
528: PetscObjectStateIncrease((PetscObject)r);
529: return(0);
530: }
534: /*@
535: PetscRandomGetValueReal - Generates a purely real random number. Call this after first calling
536: PetscRandomCreate().
538: Not Collective
540: Intput Parameter:
541: . r - the random number generator context
543: Output Parameter:
544: . val - the value
546: Level: intermediate
548: Notes:
549: Use VecSetRandom() to set the elements of a vector to random numbers.
551: Example of Usage:
552: .vb
553: PetscRandomCreate(PETSC_COMM_WORLD,&r);
554: PetscRandomGetValueReal(r,&value1);
555: PetscRandomGetValueReal(r,&value2);
556: PetscRandomGetValueReal(r,&value3);
557: PetscRandomDestroy(r);
558: .ve
560: Concepts: random numbers^getting
562: .seealso: PetscRandomCreate(), PetscRandomDestroy(), VecSetRandom(), PetscRandomGetValue()
563: @*/
564: PetscErrorCode PetscRandomGetValueReal(PetscRandom r,PetscReal *val)
565: {
573: (*r->ops->getvaluereal)(r,val);
574: PetscObjectStateIncrease((PetscObject)r);
575: return(0);
576: }