15c6c1daeSBarry Smith 25c6c1daeSBarry Smith /* 35c6c1daeSBarry Smith This file contains routines for interfacing to random number generators. 45c6c1daeSBarry Smith This provides more than just an interface to some system random number 55c6c1daeSBarry Smith generator: 65c6c1daeSBarry Smith 75c6c1daeSBarry Smith Numbers can be shuffled for use as random tuples 85c6c1daeSBarry Smith 95c6c1daeSBarry Smith Multiple random number generators may be used 105c6c1daeSBarry Smith 115c6c1daeSBarry Smith We are still not sure what interface we want here. There should be 125c6c1daeSBarry Smith one to reinitialize and set the seed. 135c6c1daeSBarry Smith */ 145c6c1daeSBarry Smith 155c6c1daeSBarry Smith #include <../src/sys/classes/random/randomimpl.h> /*I "petscsys.h" I*/ 165c6c1daeSBarry Smith 175c6c1daeSBarry Smith /*@ 185c6c1daeSBarry Smith PetscRandomGetValue - Generates a random number. Call this after first calling 195c6c1daeSBarry Smith PetscRandomCreate(). 205c6c1daeSBarry Smith 215c6c1daeSBarry Smith Not Collective 225c6c1daeSBarry Smith 235c6c1daeSBarry Smith Intput Parameter: 245c6c1daeSBarry Smith . r - the random number generator context 255c6c1daeSBarry Smith 265c6c1daeSBarry Smith Output Parameter: 275c6c1daeSBarry Smith . val - the value 285c6c1daeSBarry Smith 295c6c1daeSBarry Smith Level: intermediate 305c6c1daeSBarry Smith 315c6c1daeSBarry Smith Notes: 325c6c1daeSBarry Smith Use VecSetRandom() to set the elements of a vector to random numbers. 335c6c1daeSBarry Smith 345c6c1daeSBarry Smith When PETSc is compiled for complex numbers this returns a complex number with random real and complex parts. 355c6c1daeSBarry Smith Use PetscGetValueReal() to get a random real number. 365c6c1daeSBarry Smith 375c6c1daeSBarry Smith To get a complex number with only a random real part, first call PetscRandomSetInterval() with a equal 385c6c1daeSBarry Smith low and high imaginary part. Similarly to get a complex number with only a random imaginary part call 395c6c1daeSBarry Smith PetscRandomSetInterval() with a equal low and high real part. 405c6c1daeSBarry Smith 415c6c1daeSBarry Smith Example of Usage: 425c6c1daeSBarry Smith .vb 435c6c1daeSBarry Smith PetscRandomCreate(PETSC_COMM_WORLD,&r); 445c6c1daeSBarry Smith PetscRandomGetValue(r,&value1); 455c6c1daeSBarry Smith PetscRandomGetValue(r,&value2); 465c6c1daeSBarry Smith PetscRandomGetValue(r,&value3); 475c6c1daeSBarry Smith PetscRandomDestroy(&r); 485c6c1daeSBarry Smith .ve 495c6c1daeSBarry Smith 505c6c1daeSBarry Smith .seealso: PetscRandomCreate(), PetscRandomDestroy(), VecSetRandom(), PetscRandomGetValueReal() 515c6c1daeSBarry Smith @*/ 525c6c1daeSBarry Smith PetscErrorCode PetscRandomGetValue(PetscRandom r,PetscScalar *val) 535c6c1daeSBarry Smith { 545c6c1daeSBarry Smith PetscErrorCode ierr; 555c6c1daeSBarry Smith 565c6c1daeSBarry Smith PetscFunctionBegin; 575c6c1daeSBarry Smith PetscValidHeaderSpecific(r,PETSC_RANDOM_CLASSID,1); 585c6c1daeSBarry Smith PetscValidType(r,1); 59*808ba619SStefano Zampini if (!r->ops->getvalue) { 60*808ba619SStefano Zampini if (!r->ops->getvalues) SETERRQ1(PetscObjectComm((PetscObject)r),PETSC_ERR_SUP,"Random type %s cannot generate PetscScalar",((PetscObject)r)->type_name); 61*808ba619SStefano Zampini ierr = (*r->ops->getvalues)(r,1,val);CHKERRQ(ierr); 62*808ba619SStefano Zampini } else { 635c6c1daeSBarry Smith ierr = (*r->ops->getvalue)(r,val);CHKERRQ(ierr); 64*808ba619SStefano Zampini } 655c6c1daeSBarry Smith ierr = PetscObjectStateIncrease((PetscObject)r);CHKERRQ(ierr); 665c6c1daeSBarry Smith PetscFunctionReturn(0); 675c6c1daeSBarry Smith } 685c6c1daeSBarry Smith 695c6c1daeSBarry Smith /*@ 70*808ba619SStefano Zampini PetscRandomGetValueReal - Generates a real random number. Call this after first calling 715c6c1daeSBarry Smith PetscRandomCreate(). 725c6c1daeSBarry Smith 735c6c1daeSBarry Smith Not Collective 745c6c1daeSBarry Smith 755c6c1daeSBarry Smith Intput Parameter: 765c6c1daeSBarry Smith . r - the random number generator context 775c6c1daeSBarry Smith 785c6c1daeSBarry Smith Output Parameter: 795c6c1daeSBarry Smith . val - the value 805c6c1daeSBarry Smith 815c6c1daeSBarry Smith Level: intermediate 825c6c1daeSBarry Smith 835c6c1daeSBarry Smith Notes: 845c6c1daeSBarry Smith Use VecSetRandom() to set the elements of a vector to random numbers. 855c6c1daeSBarry Smith 865c6c1daeSBarry Smith Example of Usage: 875c6c1daeSBarry Smith .vb 885c6c1daeSBarry Smith PetscRandomCreate(PETSC_COMM_WORLD,&r); 895c6c1daeSBarry Smith PetscRandomGetValueReal(r,&value1); 905c6c1daeSBarry Smith PetscRandomGetValueReal(r,&value2); 915c6c1daeSBarry Smith PetscRandomGetValueReal(r,&value3); 925c6c1daeSBarry Smith PetscRandomDestroy(&r); 935c6c1daeSBarry Smith .ve 945c6c1daeSBarry Smith 955c6c1daeSBarry Smith .seealso: PetscRandomCreate(), PetscRandomDestroy(), VecSetRandom(), PetscRandomGetValue() 965c6c1daeSBarry Smith @*/ 975c6c1daeSBarry Smith PetscErrorCode PetscRandomGetValueReal(PetscRandom r,PetscReal *val) 985c6c1daeSBarry Smith { 995c6c1daeSBarry Smith PetscErrorCode ierr; 1005c6c1daeSBarry Smith 1015c6c1daeSBarry Smith PetscFunctionBegin; 1025c6c1daeSBarry Smith PetscValidHeaderSpecific(r,PETSC_RANDOM_CLASSID,1); 1035c6c1daeSBarry Smith PetscValidType(r,1); 104*808ba619SStefano Zampini if (!r->ops->getvaluereal) { 105*808ba619SStefano Zampini if (!r->ops->getvaluesreal) SETERRQ1(PetscObjectComm((PetscObject)r),PETSC_ERR_SUP,"Random type %s cannot generate PetscReal",((PetscObject)r)->type_name); 106*808ba619SStefano Zampini ierr = (*r->ops->getvaluesreal)(r,1,val);CHKERRQ(ierr); 107*808ba619SStefano Zampini } else { 1085c6c1daeSBarry Smith ierr = (*r->ops->getvaluereal)(r,val);CHKERRQ(ierr); 109*808ba619SStefano Zampini } 110*808ba619SStefano Zampini ierr = PetscObjectStateIncrease((PetscObject)r);CHKERRQ(ierr); 111*808ba619SStefano Zampini PetscFunctionReturn(0); 112*808ba619SStefano Zampini } 113*808ba619SStefano Zampini 114*808ba619SStefano Zampini /*@ 115*808ba619SStefano Zampini PetscRandomGetValues - Generates a sequence of random numbers. Call this after first calling 116*808ba619SStefano Zampini PetscRandomCreate(). 117*808ba619SStefano Zampini 118*808ba619SStefano Zampini Not Collective 119*808ba619SStefano Zampini 120*808ba619SStefano Zampini Intput Parameter: 121*808ba619SStefano Zampini + r - the random number generator context 122*808ba619SStefano Zampini - n - number of random numbers to generate 123*808ba619SStefano Zampini 124*808ba619SStefano Zampini Output Parameter: 125*808ba619SStefano Zampini . val - the array to hold the values 126*808ba619SStefano Zampini 127*808ba619SStefano Zampini Level: intermediate 128*808ba619SStefano Zampini 129*808ba619SStefano Zampini Notes: 130*808ba619SStefano Zampini Use VecSetRandom() to set the elements of a vector to random numbers. 131*808ba619SStefano Zampini 132*808ba619SStefano Zampini When PETSc is compiled for complex numbers this returns an array of complex numbers with random real and complex parts. 133*808ba619SStefano Zampini Use PetscGetValuesReal() to get an array of random real numbers. 134*808ba619SStefano Zampini 135*808ba619SStefano Zampini .seealso: PetscRandomCreate(), PetscRandomDestroy(), VecSetRandom(), PetscRandomGetValue() 136*808ba619SStefano Zampini @*/ 137*808ba619SStefano Zampini PetscErrorCode PetscRandomGetValues(PetscRandom r, PetscInt n, PetscScalar *val) 138*808ba619SStefano Zampini { 139*808ba619SStefano Zampini PetscErrorCode ierr; 140*808ba619SStefano Zampini 141*808ba619SStefano Zampini PetscFunctionBegin; 142*808ba619SStefano Zampini PetscValidHeaderSpecific(r,PETSC_RANDOM_CLASSID,1); 143*808ba619SStefano Zampini PetscValidType(r,1); 144*808ba619SStefano Zampini if (!r->ops->getvalues) { 145*808ba619SStefano Zampini PetscInt i; 146*808ba619SStefano Zampini if (!r->ops->getvalue) SETERRQ1(PetscObjectComm((PetscObject)r),PETSC_ERR_SUP,"Random type %s cannot generate PetscScalar",((PetscObject)r)->type_name); 147*808ba619SStefano Zampini for (i = 0; i < n; i++) { 148*808ba619SStefano Zampini ierr = (*r->ops->getvalue)(r,val+i);CHKERRQ(ierr); 149*808ba619SStefano Zampini } 150*808ba619SStefano Zampini } else { 151*808ba619SStefano Zampini ierr = (*r->ops->getvalues)(r,n,val);CHKERRQ(ierr); 152*808ba619SStefano Zampini } 153*808ba619SStefano Zampini ierr = PetscObjectStateIncrease((PetscObject)r);CHKERRQ(ierr); 154*808ba619SStefano Zampini PetscFunctionReturn(0); 155*808ba619SStefano Zampini } 156*808ba619SStefano Zampini 157*808ba619SStefano Zampini /*@ 158*808ba619SStefano Zampini PetscRandomGetValuesReal - Generates a sequence of real random numbers. Call this after first calling 159*808ba619SStefano Zampini PetscRandomCreate(). 160*808ba619SStefano Zampini 161*808ba619SStefano Zampini Not Collective 162*808ba619SStefano Zampini 163*808ba619SStefano Zampini Intput Parameter: 164*808ba619SStefano Zampini + r - the random number generator context 165*808ba619SStefano Zampini - n - number of random numbers to generate 166*808ba619SStefano Zampini 167*808ba619SStefano Zampini Output Parameter: 168*808ba619SStefano Zampini . val - the array to hold the values 169*808ba619SStefano Zampini 170*808ba619SStefano Zampini Level: intermediate 171*808ba619SStefano Zampini 172*808ba619SStefano Zampini Notes: 173*808ba619SStefano Zampini Use VecSetRandom() to set the elements of a vector to random numbers. 174*808ba619SStefano Zampini 175*808ba619SStefano Zampini .seealso: PetscRandomCreate(), PetscRandomDestroy(), VecSetRandom(), PetscRandomGetValues() 176*808ba619SStefano Zampini @*/ 177*808ba619SStefano Zampini PetscErrorCode PetscRandomGetValuesReal(PetscRandom r, PetscInt n, PetscReal *val) 178*808ba619SStefano Zampini { 179*808ba619SStefano Zampini PetscErrorCode ierr; 180*808ba619SStefano Zampini 181*808ba619SStefano Zampini PetscFunctionBegin; 182*808ba619SStefano Zampini PetscValidHeaderSpecific(r,PETSC_RANDOM_CLASSID,1); 183*808ba619SStefano Zampini PetscValidType(r,1); 184*808ba619SStefano Zampini if (!r->ops->getvaluesreal) { 185*808ba619SStefano Zampini PetscInt i; 186*808ba619SStefano Zampini if (!r->ops->getvaluereal) SETERRQ1(PetscObjectComm((PetscObject)r),PETSC_ERR_SUP,"Random type %s cannot generate PetscReal",((PetscObject)r)->type_name); 187*808ba619SStefano Zampini for (i = 0; i < n; i++) { 188*808ba619SStefano Zampini ierr = (*r->ops->getvaluereal)(r,val+i);CHKERRQ(ierr); 189*808ba619SStefano Zampini } 190*808ba619SStefano Zampini } else { 191*808ba619SStefano Zampini ierr = (*r->ops->getvaluesreal)(r,n,val);CHKERRQ(ierr); 192*808ba619SStefano Zampini } 1935c6c1daeSBarry Smith ierr = PetscObjectStateIncrease((PetscObject)r);CHKERRQ(ierr); 1945c6c1daeSBarry Smith PetscFunctionReturn(0); 1955c6c1daeSBarry Smith } 1965c6c1daeSBarry Smith 1975c6c1daeSBarry Smith /*@ 1985c6c1daeSBarry Smith PetscRandomGetInterval - Gets the interval over which the random numbers 1995c6c1daeSBarry Smith will be randomly distributed. By default, this interval is [0,1). 2005c6c1daeSBarry Smith 2015c6c1daeSBarry Smith Not collective 2025c6c1daeSBarry Smith 2035c6c1daeSBarry Smith Input Parameters: 2045c6c1daeSBarry Smith . r - the random number generator context 2055c6c1daeSBarry Smith 2065c6c1daeSBarry Smith Output Parameters: 2075c6c1daeSBarry Smith + low - The lower bound of the interval 2085c6c1daeSBarry Smith - high - The upper bound of the interval 2095c6c1daeSBarry Smith 2105c6c1daeSBarry Smith Level: intermediate 2115c6c1daeSBarry Smith 2125c6c1daeSBarry Smith .seealso: PetscRandomCreate(), PetscRandomSetInterval() 2135c6c1daeSBarry Smith @*/ 2145c6c1daeSBarry Smith PetscErrorCode PetscRandomGetInterval(PetscRandom r,PetscScalar *low,PetscScalar *high) 2155c6c1daeSBarry Smith { 2165c6c1daeSBarry Smith PetscFunctionBegin; 2175c6c1daeSBarry Smith PetscValidHeaderSpecific(r,PETSC_RANDOM_CLASSID,1); 2185c6c1daeSBarry Smith if (low) { 2195c6c1daeSBarry Smith PetscValidScalarPointer(low,2); 2205c6c1daeSBarry Smith *low = r->low; 2215c6c1daeSBarry Smith } 2225c6c1daeSBarry Smith if (high) { 2235c6c1daeSBarry Smith PetscValidScalarPointer(high,3); 2245c6c1daeSBarry Smith *high = r->low+r->width; 2255c6c1daeSBarry Smith } 2265c6c1daeSBarry Smith PetscFunctionReturn(0); 2275c6c1daeSBarry Smith } 2285c6c1daeSBarry Smith 2295c6c1daeSBarry Smith /*@ 2305c6c1daeSBarry Smith PetscRandomSetInterval - Sets the interval over which the random numbers 2315c6c1daeSBarry Smith will be randomly distributed. By default, this interval is [0,1). 2325c6c1daeSBarry Smith 2335c6c1daeSBarry Smith Not collective 2345c6c1daeSBarry Smith 2355c6c1daeSBarry Smith Input Parameters: 2365c6c1daeSBarry Smith + r - the random number generator context 2375c6c1daeSBarry Smith . low - The lower bound of the interval 2385c6c1daeSBarry Smith - high - The upper bound of the interval 2395c6c1daeSBarry Smith 2405c6c1daeSBarry Smith Level: intermediate 2415c6c1daeSBarry Smith 24295452b02SPatrick Sanan Notes: 24395452b02SPatrick Sanan 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. 2445c6c1daeSBarry Smith 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. 2455c6c1daeSBarry Smith 2465c6c1daeSBarry Smith .seealso: PetscRandomCreate(), PetscRandomGetInterval() 2475c6c1daeSBarry Smith @*/ 2485c6c1daeSBarry Smith PetscErrorCode PetscRandomSetInterval(PetscRandom r,PetscScalar low,PetscScalar high) 2495c6c1daeSBarry Smith { 2505c6c1daeSBarry Smith PetscFunctionBegin; 2515c6c1daeSBarry Smith PetscValidHeaderSpecific(r,PETSC_RANDOM_CLASSID,1); 2525c6c1daeSBarry Smith #if defined(PETSC_USE_COMPLEX) 2539f0394f4SBarry Smith if (PetscRealPart(low) > PetscRealPart(high)) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"only low <= high"); 2549f0394f4SBarry Smith if (PetscImaginaryPart(low) > PetscImaginaryPart(high)) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"only low <= high"); 2555c6c1daeSBarry Smith #else 2569f0394f4SBarry Smith if (low >= high) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"only low <= high: Instead %g %g",(double)low,(double)high); 2575c6c1daeSBarry Smith #endif 2585c6c1daeSBarry Smith r->low = low; 2595c6c1daeSBarry Smith r->width = high-low; 2605c6c1daeSBarry Smith r->iset = PETSC_TRUE; 2615c6c1daeSBarry Smith PetscFunctionReturn(0); 2625c6c1daeSBarry Smith } 263