推 gdgy:rejection 和參數轉換有何不同啊？ 05/20 23:28
推 iambbegg:方法不同 05/21 00:02
※ 編輯: iambbegg 來自: 59.104.239.175 (05/21 00:18)
→ Procyon:謝謝 05/21 08:02
※ 引述《Procyon (...)》之銘言：
: 如題
: 怎麼用C++跑出1000個為normal distribution的數字阿??
: 謝謝囉~
之前作業的一部份，這是利用rejection method產生Z~N(0,1)
若要產生Y~N(μ,σ^2)，則Y=μ+σZ (pf by mgf)
有興趣的話參考Ross的A course in simulation
-------------------------------------
//DevC++ 4.9(g++ 3.4.2)
Source code
#include<iostream>
#include<cmath>
#include<ctime>
using namespace std;
double uniform();
double rejection();
const int MAX=1000000;
int main()
{
srand(time(NULL));
double a=0;
for(int i=0;i<MAX;i++)a+=rejection();
cout<<"bar_X="<<a/MAX<<endl; //sample mean
system("pause");
return 0;
}
double rejection()//generate N(0,1) by rejection method
{
double Y;//exp(1)
double U;//uniform(0,1)
do{
Y=(-1)*log(uniform());
U=uniform();
}while(U>exp((-1)*pow(Y-1,2.)/2));
return uniform()<=0.5?Y:-Y; //distribute to Z or -Z
}
double uniform()
{
int k=(int)pow(2.,15.);
return rand()%k/(double)k;
}
------------------------------
Output
bar_X=0.00025193
--
--
※ 發信站: 批踢踢實業坊(ptt.cc)
◆ From: 59.104.239.175