Random Normal Distribution

# İlk defa Avadis Hoca'nın dersinde şu Central Limit Theorem'le karşılaşmıştım. Bize Excel'de 1000 tane uniform dağılmış random sayıyı 50'şer 50'şer toplatmış, nihayetinde elde ettiğimiz sayıların histogramını çizmemizi istemişti. Bir de ne görelim: Gaussian dağılım çıkmasın mı!
# Geçen gün Taylan Hoca SystemC öğrenmem için bir normal dağılıma uygun olarak sayı üreten bir de bu sayıların histogramını hesaplayan iki modül yazmamı istemişti. Gaussian üretecini 50'şer sayı üretip onları toplayarak mı yapacağımı sordum. Taylan Hoca da Koonin ve Meredith'in Computational Physics kitabından bir algoritma gösterdi. C kodu olarak şöyle:

double randNorm(double mean, double sigma) {
static int iset = 0;
static double gset;
double fac,rsq,v1,v2;

if (iset == 0) {
do {
v1 = randUnif(-1, 1);
v2 = randUnif(-1, 1);
rsq = v1*v1 + v2*v2;
} while (rsq >= 1.0 || rsq == 0.0);
fac = sqrt(-2.0*log(rsq)/rsq);
gset = v1*fac;
iset = 1;
return mean + v2*fac*sigma;
} else {
iset=0;
return mean + gset*sigma;
}
}

# Bunu incelemek lazım. Neden böyle. Kitabı bulup bakayım. Hımm...
# (x1,x2) civarına denk gelecek nokta sayısı: exp(-(x1^2+x2^2)/2)dx1 dx2
# polar koordinatlara geçersek: r=sqrt(x1^2+x2^2), theta=arctan(x2/x1)
# ve dağılım: exp(-r^2/2)r dr dtheta, ve -r^2/2=u dersek
# dağılım: exp(-u)dudtheta, ki u 0 ile sonsuz aradında, theta'da 0 ile 2pi arasında değişiyor.
# Bu bağlamda, x1=rcos(theta)=(2u)^(1/2)cos(theta), x2=(2u)^(1/2)sin(theta) normal dağılacak.
# Kitaptaki BASIC kodu:

TWOU = -2*LOG(1-RND)
RADIUS=SQR(TWOU)
THETA=2*3.14159*RND
GAUSS1=RADIUS*COS(THETA)
GAUSS2=RADIUS*SIN(THETA)
RETURN

# İstenen ortalama ve sigma'da bir dağılım elde etmek için de üretilen sayıyı sigma ile çarpıp xbar eklemek yetiyormuş. Bu daha anlaşılır geldi şimdi bana.