ページ 11

rand関数 課題

Posted: 2011年7月01日(金) 01:23
by SEED
問題
1cmの正方形の中に半径1cmの円があり、それを1/4にする。
この図形に球をたくさん落とす。
N回落としたとき、扇形の内側に落ちたときの数をnとする。
球の座標を(x,y)として、xとyをrand関数で取り出し、n/Nとしπの近似値を求めよ。

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

int main(void){

float x,y,pi;
int N,n=0;
int i;
N=10;
for(i=0;i<N;i++){

srand((unsigned) time(NULL));
x=rand()%100000/100000.0;
y=rand()%100000/100000.0;

if(x*x+y*y<=1)
n++;
}

pi=4.0*n/N;
printf("π=%f\n",pi);

return 0;
}

ここまで作ったのですが、Nの値を何にしても結果が4になってしまいます。
よければ解決法を教えて貰えるとありがたいです。

Re: rand関数 課題

Posted: 2011年7月01日(金) 01:30
by h2so5
おそらくsrandをループ内で毎回呼んでいるのが原因です。

time関数は秒単位の時間を返しますので、
ループに1秒もかからない場合乱数が毎回同じ数値でリセットされてしまい、
rand()の返す値が変わりません。

Re: rand関数 課題

Posted: 2011年7月01日(金) 04:01
by amane
これはモンテカルロ法ですね。
ソース自体に問題はないと思いますが……

少し関係ありませんが、まず問題のイメージが出来ていますか?
>1cmの正方形の中に半径1cmの円があり、それを1/4にする。
もしかして、半径1cm円の1/4の扇形が、1cmの正方形の中に入っているのではないですか?
1cmの正方形の中に半径1cmの円は入りませんね。
正方形の左下を中心に、半径1cmの1/4円があるとイメージしてみてください。

「Nの値を何にしても」というのは、どの程度の範囲ででしょうか?
100、1000、10000などの大きい数にしてみましたか?
もしかしたら、Nの値が小さいことが原因かもしれません。
正方形に球を落として、扇内に入った数をnとするので、
落とす数があまりに少ないと、ほとんどが扇内に落ちることがイメージできますよね。

Re: rand関数 課題

Posted: 2011年7月01日(金) 09:25
by non
SEEDさんが使っているコンパイラのRAND_MAX はいくらですか?

Re: rand関数 課題

Posted: 2011年7月01日(金) 09:28
by さかまき
>x=rand()%100000/100000.0;
>y=rand()%100000/100000.0;
実際にどんな値がx,yにセットされているかを調べてみましょう。

Re: rand関数 課題

Posted: 2011年7月01日(金) 09:35
by さかまき
すべて同じ値になるのは、srand((unsigned) time(NULL));の位置。
すべてx*x+y*y<=1になってしまうのは、
rand()%100000/100000.0 の処理とRAND_MAXとの関係からです。

Re: rand関数 課題

Posted: 2011年7月01日(金) 17:01
by GRAM
さかまき さんが書きました:>x=rand()%100000/100000.0;
>y=rand()%100000/100000.0;
実際にどんな値がx,yにセットされているかを調べてみましょう。
一応補足しておきますが、ややこしいことにおそらくVisualStudioのステップ実行などを使ってしまうと一見正常に見えてしまうことでしょう
(あくまでtime(NULL));が短時間の間にたくさん呼び出されることが問題となっているため)
また環境によってはrand関数の位置を調整しても100000もループすればπの値に収束しないこともあり得ます。
(rand関数の乱数の質の問題により)

randの質に目をつむった場合、簡単な解決策はsrandをforループの外に出してしまうことです。
3.141くらいまで収束していれば上出来じゃないかと思います。

追記:
>>nonさん
あぁ確かに引用もとをRAND_MAXについての指摘としてみた場合誤認といえますね。
自分はtimeの話だとみた場合xとyの値をprintfとかで実際に出力しないとわからない、という指摘をしたつもりだったんですけどね
いやはや。

Re: rand関数 課題

Posted: 2011年7月01日(金) 20:37
by non
GRAMさん。さかまきさんが指摘していることを誤認しています。(恐らく)

Re: rand関数 課題

Posted: 2011年7月01日(金) 22:02
by box
GRAM さんが書きました: また環境によってはrand関数の位置を調整しても100000もループすればπの値に収束しないこともあり得ます。
何回読み返してもこの文の意味が分からないです。
特に「100000もループすればπの値に収束しないこともあり得ます」のところが。

Re: rand関数 課題

Posted: 2011年7月02日(土) 00:06
by ISLe
さかまきさんが指摘しているとおりrand()は0~RAND_MAXの値を返します。
そして多くの処理系でRAND_MAXは32767です。

なので
x=rand()%100000/100000.0;
y=rand()%100000/100000.0;
というコードで、座標(x,y)が円周外になることは決して無いわけですね。