モンテカルロ法の面積近似
Mが2の32乗の、線形合同法で乱数を発生し、
モンテカルロ法を使って、
領域a = {(x,y)|x≧1,y≧2,(x-1)^2 + (y-2)^2 ≦1}の面積の近似値を計算するプログラムを作成し、これを利用してn=10000とした場合の上の領域の面積、および円周率πの近似値を求めよ。
(関数を使う)
という課題が出たのですが、乱数の出し方？が変であまり近似されません(^_^;)
どなたかどこが変か教えてください(>_<)
これだと、2と4になってしまうんです(/_;)
#include<stdio.h>
#include<stdlib.h>
#include<math.h>
#define N 100
/*領域Rを含む面積T 領域R=求める面積*/
#define Nx 1
#define Ny 2
/*乱数を発生させる関数*/
unsigned long int random_g(unsigned long int x0, unsigned long int a, unsigned long int c);
int main(){
unsigned long int x0, a, c, M;
double x, y, s;
int n, count_in, count_out;
count_in = 0;
count_out = 0;
printf("x0: ", x0 ); scanf("%d", &x0);
a = 61;
c = 49;
for(n = 0; n < 10000; n++){
/*乱数を発生*/
x = (random_g(x0, 4*a*n+1, c) / (double)4294967296)*Nx + 1;
y = (random_g(x0, 4*a*n+1, c) / (double)4294967296)*Ny + 1;
if(x >= 1 && y >= 1 && pow(x-1, 2) + pow(y-2, 2) <= 1)
count_in++;
else
count_out++;
}
s = (double)count_in/(count_in + count_out)*Nx*Ny;
printf("s = %f\n", s);
printf("pi = %f\n", 2*s);
return 0;
}
unsigned long int random_g(unsigned long int x0, unsigned long int a, unsigned long int c){
int i;
x0 = (a*x0 + c);
/*printf("%u\n", x0);*/
return x0;
}
お礼
わかりました。 詳しい解説ありがとうございます。