1 /* 2 This file contains routines for interfacing to random number generators. 3 This provides more than just an interface to some system random number 4 generator: 5 6 Numbers can be shuffled for use as random tuples 7 8 Multiple random number generators may be used 9 10 We are still not sure what interface we want here. There should be 11 one to reinitialize and set the seed. 12 */ 13 14 #include <petsc/private/randomimpl.h> /*I "petscsys.h" I*/ 15 #include <petscviewer.h> 16 17 /* Logging support */ 18 PetscClassId PETSC_RANDOM_CLASSID; 19 20 /*@ 21 PetscRandomDestroy - Destroys a `PetscRandom` object that was created by `PetscRandomCreate()`. 22 23 Collective 24 25 Input Parameter: 26 . r - the random number generator object 27 28 Level: intermediate 29 30 .seealso: `PetscRandom`, `PetscRandomGetValue()`, `PetscRandomCreate()`, `VecSetRandom()` 31 @*/ 32 PetscErrorCode PetscRandomDestroy(PetscRandom *r) 33 { 34 PetscFunctionBegin; 35 if (!*r) PetscFunctionReturn(PETSC_SUCCESS); 36 PetscValidHeaderSpecific(*r, PETSC_RANDOM_CLASSID, 1); 37 if (--((PetscObject)*r)->refct > 0) { 38 *r = NULL; 39 PetscFunctionReturn(PETSC_SUCCESS); 40 } 41 PetscTryTypeMethod(*r, destroy); 42 PetscCall(PetscHeaderDestroy(r)); 43 PetscFunctionReturn(PETSC_SUCCESS); 44 } 45 46 /*@C 47 PetscRandomGetSeed - Gets the random seed. 48 49 Not collective 50 51 Input Parameter: 52 . r - The random number generator context 53 54 Output Parameter: 55 . seed - The random seed 56 57 Level: intermediate 58 59 .seealso: `PetscRandom`, `PetscRandomCreate()`, `PetscRandomSetSeed()`, `PetscRandomSeed()` 60 @*/ 61 PetscErrorCode PetscRandomGetSeed(PetscRandom r, PetscInt64 *seed) 62 { 63 PetscFunctionBegin; 64 PetscValidHeaderSpecific(r, PETSC_RANDOM_CLASSID, 1); 65 if (seed) { 66 PetscAssertPointer(seed, 2); 67 *seed = (PetscInt64)r->seed; 68 } 69 PetscFunctionReturn(PETSC_SUCCESS); 70 } 71 72 /*@C 73 PetscRandomSetSeed - Sets the random seed. You MUST call `PetscRandomSeed()` after this call to have the new seed used. 74 75 Not collective 76 77 Input Parameters: 78 + r - The random number generator context 79 - seed - The random seed 80 81 Level: intermediate 82 83 Example Usage: 84 .vb 85 PetscRandomSetSeed(r,a positive integer); 86 PetscRandomSeed(r); 87 PetscRandomGetValue() will now start with the new seed. 88 89 PetscRandomSeed(r) without a call to PetscRandomSetSeed() re-initializes 90 the seed. The random numbers generated will be the same as before. 91 .ve 92 93 .seealso: `PetscRandom`, `PetscRandomCreate()`, `PetscRandomGetSeed()`, `PetscRandomSeed()` 94 @*/ 95 PetscErrorCode PetscRandomSetSeed(PetscRandom r, PetscInt64 seed) 96 { 97 PetscFunctionBegin; 98 PetscValidHeaderSpecific(r, PETSC_RANDOM_CLASSID, 1); 99 r->seed = (unsigned long)seed; 100 PetscCall(PetscInfo(NULL, "Setting seed to %d\n", (int)seed)); 101 PetscFunctionReturn(PETSC_SUCCESS); 102 } 103 104 /* ------------------------------------------------------------------- */ 105 /* 106 PetscRandomSetTypeFromOptions_Private - Sets the type of random generator from user options. Defaults to type PETSCRAND48 or PETSCRAND. 107 108 Collective 109 110 Input Parameter: 111 . rnd - The random number generator context 112 113 Level: intermediate 114 115 .seealso: `PetscRandomSetFromOptions()`, `PetscRandomSetType()` 116 */ 117 static PetscErrorCode PetscRandomSetTypeFromOptions_Private(PetscRandom rnd, PetscOptionItems PetscOptionsObject) 118 { 119 PetscBool opt; 120 const char *defaultType; 121 char typeName[256]; 122 123 PetscFunctionBegin; 124 if (((PetscObject)rnd)->type_name) { 125 defaultType = ((PetscObject)rnd)->type_name; 126 } else { 127 defaultType = PETSCRANDER48; 128 } 129 130 PetscCall(PetscRandomRegisterAll()); 131 PetscCall(PetscOptionsFList("-random_type", "PetscRandom type", "PetscRandomSetType", PetscRandomList, defaultType, typeName, 256, &opt)); 132 if (opt) { 133 PetscCall(PetscRandomSetType(rnd, typeName)); 134 } else { 135 PetscCall(PetscRandomSetType(rnd, defaultType)); 136 } 137 PetscFunctionReturn(PETSC_SUCCESS); 138 } 139 140 /*@ 141 PetscRandomSetFromOptions - Configures the random number generator from the options database. 142 143 Collective 144 145 Input Parameter: 146 . rnd - The random number generator context 147 148 Options Database Keys: 149 + -random_seed <integer> - provide a seed to the random number generator 150 - -random_no_imaginary_part - makes the imaginary part of the random number zero, this is useful when you want the 151 same code to produce the same result when run with real numbers or complex numbers for regression testing purposes 152 153 Level: beginner 154 155 Note: 156 Must be called after `PetscRandomCreate()` but before the rnd is used. 157 158 .seealso: `PetscRandom`, `PetscRandomCreate()`, `PetscRandomSetType()` 159 @*/ 160 PetscErrorCode PetscRandomSetFromOptions(PetscRandom rnd) 161 { 162 PetscBool set, noimaginary = PETSC_FALSE; 163 PetscInt seed; 164 165 PetscFunctionBegin; 166 PetscValidHeaderSpecific(rnd, PETSC_RANDOM_CLASSID, 1); 167 168 PetscObjectOptionsBegin((PetscObject)rnd); 169 170 /* Handle PetscRandom type options */ 171 PetscCall(PetscRandomSetTypeFromOptions_Private(rnd, PetscOptionsObject)); 172 173 /* Handle specific random generator's options */ 174 PetscTryTypeMethod(rnd, setfromoptions, PetscOptionsObject); 175 PetscCall(PetscOptionsInt("-random_seed", "Seed to use to generate random numbers", "PetscRandomSetSeed", 0, &seed, &set)); 176 if (set) { 177 PetscCall(PetscRandomSetSeed(rnd, (unsigned long int)seed)); 178 PetscCall(PetscRandomSeed(rnd)); 179 } 180 PetscCall(PetscOptionsBool("-random_no_imaginary_part", "The imaginary part of the random number will be zero", "PetscRandomSetInterval", noimaginary, &noimaginary, &set)); 181 #if defined(PETSC_HAVE_COMPLEX) 182 if (set) { 183 if (noimaginary) { 184 PetscScalar low, high; 185 PetscCall(PetscRandomGetInterval(rnd, &low, &high)); 186 low = low - PetscImaginaryPart(low); 187 high = high - PetscImaginaryPart(high); 188 PetscCall(PetscRandomSetInterval(rnd, low, high)); 189 } 190 } 191 #endif 192 PetscOptionsEnd(); 193 PetscCall(PetscRandomViewFromOptions(rnd, NULL, "-random_view")); 194 PetscFunctionReturn(PETSC_SUCCESS); 195 } 196 197 /*@ 198 PetscRandomSetOptionsPrefix - Sets the prefix used for searching for all 199 `PetscRandom` options in the database. 200 201 Logically Collective 202 203 Input Parameters: 204 + r - the random number generator context 205 - prefix - the prefix to prepend to all option names 206 207 Level: advanced 208 209 Note: 210 A hyphen (-) must NOT be given at the beginning of the prefix name. 211 The first character of all runtime options is AUTOMATICALLY the hyphen. 212 213 .seealso: `PetscRandom`, `PetscRandomSetFromOptions()` 214 @*/ 215 PetscErrorCode PetscRandomSetOptionsPrefix(PetscRandom r, const char prefix[]) 216 { 217 PetscFunctionBegin; 218 PetscValidHeaderSpecific(r, PETSC_RANDOM_CLASSID, 1); 219 PetscCall(PetscObjectSetOptionsPrefix((PetscObject)r, prefix)); 220 PetscFunctionReturn(PETSC_SUCCESS); 221 } 222 223 #if defined(PETSC_HAVE_SAWS) 224 #include <petscviewersaws.h> 225 #endif 226 227 /*@ 228 PetscRandomViewFromOptions - View a `PetscRandom` object based on the options database 229 230 Collective 231 232 Input Parameters: 233 + A - the random number generator context 234 . obj - Optional object 235 - name - command line option 236 237 Level: intermediate 238 239 .seealso: `PetscRandom`, `PetscRandomView`, `PetscObjectViewFromOptions()`, `PetscRandomCreate()` 240 @*/ 241 PetscErrorCode PetscRandomViewFromOptions(PetscRandom A, PetscObject obj, const char name[]) 242 { 243 PetscFunctionBegin; 244 PetscValidHeaderSpecific(A, PETSC_RANDOM_CLASSID, 1); 245 PetscCall(PetscObjectViewFromOptions((PetscObject)A, obj, name)); 246 PetscFunctionReturn(PETSC_SUCCESS); 247 } 248 249 /*@ 250 PetscRandomView - Views a random number generator object. 251 252 Collective 253 254 Input Parameters: 255 + rnd - The random number generator context 256 - viewer - an optional visualization context 257 258 Level: beginner 259 260 Note: 261 The available visualization contexts include 262 + `PETSC_VIEWER_STDOUT_SELF` - standard output (default) 263 - `PETSC_VIEWER_STDOUT_WORLD` - synchronized standard 264 output where only the first processor opens 265 the file. All other processors send their 266 data to the first processor to print. 267 268 .seealso: `PetscRandom`, `PetscRealView()`, `PetscScalarView()`, `PetscIntView()` 269 @*/ 270 PetscErrorCode PetscRandomView(PetscRandom rnd, PetscViewer viewer) 271 { 272 PetscBool iascii; 273 #if defined(PETSC_HAVE_SAWS) 274 PetscBool issaws; 275 #endif 276 277 PetscFunctionBegin; 278 PetscValidHeaderSpecific(rnd, PETSC_RANDOM_CLASSID, 1); 279 PetscValidType(rnd, 1); 280 if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)rnd), &viewer)); 281 PetscValidHeaderSpecific(viewer, PETSC_VIEWER_CLASSID, 2); 282 PetscCheckSameComm(rnd, 1, viewer, 2); 283 PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii)); 284 #if defined(PETSC_HAVE_SAWS) 285 PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERSAWS, &issaws)); 286 #endif 287 if (iascii) { 288 PetscMPIInt rank; 289 PetscCall(PetscObjectPrintClassNamePrefixType((PetscObject)rnd, viewer)); 290 PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)rnd), &rank)); 291 PetscCall(PetscViewerASCIIPushSynchronized(viewer)); 292 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "[%d] Random type %s, seed %lu\n", rank, ((PetscObject)rnd)->type_name, rnd->seed)); 293 PetscCall(PetscViewerFlush(viewer)); 294 PetscCall(PetscViewerASCIIPopSynchronized(viewer)); 295 #if defined(PETSC_HAVE_SAWS) 296 } else if (issaws) { 297 PetscMPIInt rank; 298 const char *name; 299 300 PetscCall(PetscObjectGetName((PetscObject)rnd, &name)); 301 PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank)); 302 if (!((PetscObject)rnd)->amsmem && rank == 0) { 303 char dir[1024]; 304 305 PetscCall(PetscObjectViewSAWs((PetscObject)rnd, viewer)); 306 PetscCall(PetscSNPrintf(dir, 1024, "/PETSc/Objects/%s/Low", name)); 307 PetscCallSAWs(SAWs_Register, (dir, &rnd->low, 1, SAWs_READ, SAWs_DOUBLE)); 308 } 309 #endif 310 } 311 PetscFunctionReturn(PETSC_SUCCESS); 312 } 313 314 /*@ 315 PetscRandomCreate - Creates an object for generating random numbers, 316 and initializes the random-number generator. 317 318 Collective 319 320 Input Parameter: 321 . comm - MPI communicator 322 323 Output Parameter: 324 . r - the random number generator object 325 326 Level: intermediate 327 328 Notes: 329 The random type has to be set by `PetscRandomSetType()`. 330 331 This is only a primitive "parallel" random number generator, it should NOT 332 be used for sophisticated parallel Monte Carlo methods since it will very likely 333 not have the correct statistics across processors. You can provide your own 334 parallel generator using `PetscRandomRegister()`; 335 336 If you create a `PetscRandom()` using `PETSC_COMM_SELF` on several processors then 337 the SAME random numbers will be generated on all those processors. Use `PETSC_COMM_WORLD` 338 or the appropriate parallel communicator to eliminate this issue. 339 340 Use `VecSetRandom()` to set the elements of a vector to random numbers. 341 342 Example of Usage: 343 .vb 344 PetscRandomCreate(PETSC_COMM_SELF,&r); 345 PetscRandomSetType(r,PETSCRAND48); 346 PetscRandomGetValue(r,&value1); 347 PetscRandomGetValueReal(r,&value2); 348 PetscRandomDestroy(&r); 349 .ve 350 351 .seealso: `PetscRandomSetType()`, `PetscRandomGetValue()`, `PetscRandomGetValueReal()`, `PetscRandomSetInterval()`, 352 `PetscRandomDestroy()`, `VecSetRandom()`, `PetscRandomType`, `PetscRandom` 353 @*/ 354 PetscErrorCode PetscRandomCreate(MPI_Comm comm, PetscRandom *r) 355 { 356 PetscRandom rr; 357 PetscMPIInt rank; 358 359 PetscFunctionBegin; 360 PetscAssertPointer(r, 2); 361 PetscCall(PetscRandomInitializePackage()); 362 363 PetscCall(PetscHeaderCreate(rr, PETSC_RANDOM_CLASSID, "PetscRandom", "Random number generator", "Sys", comm, PetscRandomDestroy, PetscRandomView)); 364 PetscCallMPI(MPI_Comm_rank(comm, &rank)); 365 rr->data = NULL; 366 rr->low = 0.0; 367 rr->width = 1.0; 368 rr->iset = PETSC_FALSE; 369 rr->seed = 0x12345678 + 76543 * rank; 370 PetscCall(PetscRandomSetType(rr, PETSCRANDER48)); 371 *r = rr; 372 PetscFunctionReturn(PETSC_SUCCESS); 373 } 374 375 /*@ 376 PetscRandomSeed - Seed the random number generator. 377 378 Not collective 379 380 Input Parameter: 381 . r - The random number generator context 382 383 Level: intermediate 384 385 Example Usage: 386 .vb 387 PetscRandomSetSeed(r,a positive integer); 388 PetscRandomSeed(r); 389 PetscRandomGetValue() will now start with the new seed. 390 391 PetscRandomSeed(r) without a call to PetscRandomSetSeed() re-initializes 392 the seed. The random numbers generated will be the same as before. 393 .ve 394 395 .seealso: `PetscRandomCreate()`, `PetscRandomGetSeed()`, `PetscRandomSetSeed()` 396 @*/ 397 PetscErrorCode PetscRandomSeed(PetscRandom r) 398 { 399 PetscFunctionBegin; 400 PetscValidHeaderSpecific(r, PETSC_RANDOM_CLASSID, 1); 401 PetscValidType(r, 1); 402 403 PetscUseTypeMethod(r, seed); 404 PetscCall(PetscObjectStateIncrease((PetscObject)r)); 405 PetscFunctionReturn(PETSC_SUCCESS); 406 } 407