ボックスミューラーのプログラム

ボックスミューラーのプログラムです。
アルゴリズムは下記の記事です。


ボックスミュラーのアルゴリズム - nishiru3の日記

平均、分散、個数を入力し、
発生させ、標準出力するようにしています。

ボックスミューラーで標準正規分布を作成し、線形変換より、
指定した平均と分散になるようにしています。

出力値のヒストグラムを作れば、妥当性をチェックできると思います。

C/C++で書いてます。
本当であれば、ポインタを返すようにすれば良いですが・・・。

#include <iostream>
#include <cstdlib>
#include <cmath>
using namespace std;
const int pi = 3.141592;
void boxmuller (double mu,double sigma,int imax) {
    double num1,num2,x1,x2,z1,z2;
    for (int i = 0; i < imax; i++){
	// 0~1の一様乱数生成
	num1 = (double)rand()/(double)RAND_MAX;
	num2 = (double)rand()/(double)RAND_MAX;
	// ボックスミュラー法
	x1 = sqrt(-2.0*log(num1))*cos(2*pi*num2);
	x2 = sqrt(-2.0*log(num1))*sin(2*pi*num2);
	// 線形変換
	z1  = mu + sigma*x1;
	z2  = mu + sigma*x2;
	cout << z1 << "," << z2 << endl; 
    }
}
int main(void) {
    // ボックスミューラー法による正規分布の作成
    // 平均値、分散、個数
    boxmuller(0,1,1000);	
}