C言語関係掲示板

過去ログ

No.985 モンテカルロ法で半径1のN次元の球の体積

[戻る] [ホームページ]
No.12843

モンテカルロ法で半径1のN次元の球の体積
投稿者---けいじ(2004/02/19 13:57:15)


下のプログラムを参考にして半径1のN次元の球の体積を求める方法を
教えてください。できればコメントも付けていただけたら御の字です。


/* モンテカルロ法による円の面積の計算 */

#include <stdio.h>
#include <stdlib.h>
#include <math.h>

main()
{
unsigned int seed;
int ntrial;
int i, j, nin;
double a1;
double r0, r02, r2, x, y, s;


printf("Input: Number of trials:\n");
scanf("%u", &ntrial);
printf("Input: seed\n");
scanf("%u", &seed);
/* printf("seed=%u\n",seed); */

srand(seed);

r0=1.0; /* 円の半径*/
r02=r0*r0;

nin=0;

for(i=0; i < ntrial; i++){
j=rand();
a1=j/(RAND_MAX+1.0);
x=a1;
j=rand();
a1=j/(RAND_MAX+1.0);
y=a1;
r2=x*x+y*y; /* 原点からの距離の2乗 */
if(r2 < r02) nin=nin+1; /* 円の中 */
}

s=4.0 * ((double)nin/(double)ntrial);
printf("ntrial=%u, s=%f\n", ntrial, s);

return;
}



No.12849

Re:モンテカルロ法で半径1のN次元の球の体積
投稿者---たか(2004/02/19 16:43:28)


>下のプログラムを参考にして半径1のN次元の球の体積を求める方法を
>教えてください。できればコメントも付けていただけたら御の字です。

ちょっと前に私が書いたのですがまだ過去ログ化してないようです。
その時のソースをそのまま貼り付けます。

#include <stdio.h>
#include <stdlib.h>
#include <math.h>

#define PI 3.1415926535897932384626

double montecarlo(int d, int n);
double hypersphere_volume(int n);
double fact(int n);

int main(void)
{
  int d, n;
  unsigned int seed;
  char dummy[100];
  double v, vv;
  
  printf("次元の数 D = ");
  scanf("%d", &d);
  printf("点の数 N = ");
  scanf("%d", &n);
  printf("乱数の種 SEED = ");
  scanf("%u", &seed);
  gets(dummy);
  
  srand(seed);

  v = montecarlo(d, n);
  vv = hypersphere_volume(d);

  printf("モンテカルロ法による体積 V = %fcm^%d\n計算による体積 V' = %fcm^%d\n", v, d, vv, d);
  
  return 0;
}

double montecarlo(int d, int n)
{
  int i, j, count = 0;
  double dist;
  
  for (i = 0; i < n; i++) {
    dist = 0;
    for (j = 0; j < d; j++)
      dist += pow((double)rand() / (RAND_MAX + 0.000001), 2);
    if (dist <= 1)
      count++;
  }
  return pow(2, d) * count / n;
}

double hypersphere_volume(int n)
{
  if ((n % 2) == 0)
    return pow(PI, n / 2) / fact(n / 2);
  else
    return fact((n - 1) / 2) * pow(2, n) / fact(n) * pow(PI, (n - 1) / 2);
}

double fact(int n)
{
  double res = 1;
  
  for (; n; n--)
    res *= n;
    
  return res;
}


No.12860

Re:モンテカルロ法で半径1のN次元の球の体積
投稿者---けいじ(2004/02/20 15:11:44)


>>下のプログラムを参考にして半径1のN次元の球の体積を求める方法を
>>教えてください。できればコメントも付けていただけたら御の字です。
>
>ちょっと前に私が書いたのですがまだ過去ログ化してないようです。
>その時のソースをそのまま貼り付けます。
>
><pre>#include <stdio.h>
#include <stdlib.h>
#include <math.h>

#define PI 3.1415926535897932384626

double montecarlo(int d, int n);   /*ここ*/
double hypersphere_volume(int n);  /*ここ*/
double fact(int n);          /*ここ*/

int main(void)
{
int d, n;
unsigned int seed;
char dummy[100];    /*ここ*/
double v, vv;

printf("次元の数 D = ");
scanf("%d", &d);
printf("点の数 N = ");
scanf("%d", &n);
printf("乱数の種 SEED = ");
scanf("%u", &seed);
gets(dummy);

srand(seed);

v = montecarlo(d, n);
vv = hypersphere_volume(d);

printf("モンテカルロ法による体積 V = %fcm^%d\n計算による体積 V' = %fcm^%d\n", v, d, vv, d);

return 0;
}

double montecarlo(int d, int n)
{
int i, j, count = 0;
double dist;   /*このdistのいみは?*/

for (i = 0; i < n; i++) {
dist = 0;
for (j = 0; j < d; j++)
dist += pow((double)rand() / (RAND_MAX + 0.000001), 2);
if (dist <= 1)
count++;
}
return pow(2, d) * count / n;
}

double hypersphere_volume(int n)
{
if ((n % 2) == 0)
return pow(PI, n / 2) / fact(n / 2);
else
return fact((n - 1) / 2) * pow(2, n) / fact(n) * pow(PI, (n - 1) / 2);
}

double fact(int n)  /*factの意味は?*/
{
double res = 1;

for (; n; n--)
res *= n;

return res;
}
</pre>


たかさん返信ありがとうございます。
少し理解ができないところがありますので上にコメントとして
質問を書いていますので教えていただけませんか?
おねがいします。


No.12861

Re:モンテカルロ法で半径1のN次元の球の体積
投稿者---けいじ(2004/02/20 15:42:41)


たかさんあと申し訳ないのですがこれをグラフにするために
実行結果を試行回数すべて表示するためにはどうしたらよいのでしょうか?

No.12864

Re:モンテカルロ法で半径1のN次元の球の体積
投稿者---たか(2004/02/20 16:40:37)


>たかさんあと申し訳ないのですがこれをグラフにするために
>実行結果を試行回数すべて表示するためにはどうしたらよいのでしょうか?

これなのですが、例えば試行回数が 10000回だったら、軸に10000個の
要素を取るという意味なのでしょうか。もう少し詳しくお願いします。

No.12865

Re:モンテカルロ法で半径1のN次元の球の体積
投稿者---けいじ(2004/02/20 17:09:57)


>>たかさんあと申し訳ないのですがこれをグラフにするために
>>実行結果を試行回数すべて表示するためにはどうしたらよいのでしょうか?
>
>これなのですが、例えば試行回数が 10000回だったら、軸に10000個の
>要素を取るという意味なのでしょうか。もう少し詳しくお願いします。

>グラフにするのはそういう意味です。
>あと、コマンドライン上で10000回だったらその分の値の変わりようを表示したいのです。これは、一応自分で考えてなんとか表示できるようになりました。


No.12863

Re:モンテカルロ法で半径1のN次元の球の体積
投稿者---たか(2004/02/20 16:39:31)


>double montecarlo(int d, int n);   /*ここ*/

モンテカルロ法によりd次元球をn回の試行により求める。

>double hypersphere_volume(int n);  /*ここ*/

n次元球の体積を公式により求める。(モンテカルロ法の結果との
比較用)

>double fact(int n);          /*ここ*/

fact()のプロトタイプ宣言

> char dummy[100];    /*ここ*/

入力にscanf()を使っているので、バッファに残った'\n'を”食う”
為に使っている。

> double dist;   /*このdistのいみは?*/

乱数による試行の結果が半径1を超えるかどうかの判定に使用。

>double fact(int n)  /*factの意味は?*/

階乗 n! を求める関数。

No.12866

Re:モンテカルロ法で半径1のN次元の球の体積
投稿者---けいじ(2004/02/20 17:10:38)


了解しました。
返信ありがとうございます。