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