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