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