xref: /petsc/src/sys/classes/random/interface/random.c (revision 808ba619248ec38c62492271953950a08e6d0b8c)
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