Самая простая последовательность, которую можно предложить для реализации генератора равномерного распределения:
I(j+1)=a*I(j)(mod m)
при соответствующем выборе констант. Константы были предложены Park и Miller:
a=75=16807, m=231-1=2147483647
и протестированы в исследованиях Lewis, Goodman, Miller (1969).
Прямое приложение этого метода возможно на языках ассемблера, но языки высокого уровня могут при этом зафиксировать переполнение. Для обхода этого Scharge предложил метод частичной факторизации модуля. Модуль разлагается в выражение:
m=a*q+r
Если r<q и 0<z<m-1, то при этом величины a*(z mod q) и r*[z/q] всегда лежат в интервале 0,...,m-1. Для умножения (a*z)(mod m) при этом используется алгоритм:
t = a(z mod q)-r[z/q]
если t<0, то t += m.
(a*z)(mod m)=t.
В случае констант Парка-Миллера можно использовать q=12773 и r=2836.
/* Minimal portable random generator by Park and Miller */
/* Lewis-Goodman-Miller constants */
#define IA 16807
#define IM 2147483647
#define AM (1./IM)
/* Scharge constants */
#define IQ 12773
#define IR 2836
/* Special mask to be explained below */
#define MASK 123456789
static long dummy;
/* initial seed, for all the generators here */
void Seed(long dum) {dummy=dum;}
/* returns random uniformly distributed between 0 and 1 */
float unirand0(void) {
long k;
float ans;
dummy^=MASK; /* avoid dummy==0 */
k=dummy/IQ;
if((dummy=IA*(dummy-k*IQ)-IR*k)<0) dummy+=IM;
ans=AM*dummy;
dummy^=MASK; /* restore unmasked dummy */
return(ans);
}
Использование маски связано с тем, что специфика алгоритма не позволяет устанавливать счетчик в нуль. Но, как показывает опыт, большинство пользователей счетчиков делают именно так. Маска гарантирует, что установленный счетчик не будет нулем. Если вы очень уверены в том, что человек не допустит подобной ошибки после вашего предупреждения,то можете убрать из программы все инструкции, связанные с маской.
Более продвинутая версия генератора Парка-Миллера.
В дополнение к предыдущему алгоритму, производит перетасовку по методу Bays и Durham, что позволяет уничтожить нежелательные корреляции между сериями сгенерированных случайных чисел.
/* Minimal random generator with Bays-Durham shuffle */
/* previous functions, variables and constants are described above.
Only additional follow */
#define NTAB 32
#define NWUP 8
#define NDIV (1+(IM-1)/NTAB)
#define EPS 1.2e-7
#define RNMX (1.0-EPS)
float unirand1(void) {
int j;
long k;
static long iy=0,iv[NTAB];
float temp;
/* initialize */
if(dummy<=0 || !iy) {
/* avoid negative or zero seed */
if(dummy<0) dummy=-dummy;
else if(dummy==0) dummy=1;
/* after NWUP warmups, initialize shuffle table */
for(j=NTAB+NWUP-1;j>=0;j--) {
k=dummy/IQ;
if((dummy=IA*(dummy-k*IQ)-IR*k)<0) dummy+=IM;
if(j<NTAB) iv[j]=dummy;
}
/* first specimen from the table */
iy=iv[0];
}
/* regular work: generate new number */
k=dummy/IQ;
if((dummy=IA*(dummy-k*IQ)-IR*k)<0) dummy+=IM;
/* shuffle output */
iy=iv[j=iy/NDIV];iv[j]=dummy;
/* return */
if((temp=AM*iy)>RNMX) return(RNMX);
else return(temp);
}
Заметим, что маска здесь не применяется, поскольку метод инициализации алгоритма отсекает ошибочные значения текущего счетчика.
Алгоритм Парка-Миллера с перетасовкой проходит все известные тесты, включая и те, где без перетасовки этот алгоритм дает сбой. Хорошее поведение наблюдается до значений числа вызовов счетчика порядка 108, где он начинает повторяться.
8 8 8
| 
*
Обратная связь