それ，numpy で書かない？--4--

それ，numpy で書かない？--4--

for ループが遅いのに，それを使わざるを得ない選択肢を採っていないかな?

課題：random.random() とか random.uniform() を使うものだから for ループを使わないといけなくなるというジレンマ

https://qiita.com/yamadasuzaku/items/81c9750cf96a844aec38 のモンテカルロ法によるシミュレーション部分は以下のようになっている。

    # 円と長方形の両方に含まれる点の数を数える
    overlap_count = 0
    for _ in range(num_points):
        # 一様な密度で円内の点を生成
        angle = random.uniform(0, 2 * math.pi)
        r = math.sqrt(random.uniform(0, 1)) * radius
        x = r * math.cos(angle)
        y = r * math.sin(angle)

        # この点が長方形の内部にも含まれるかチェック
        if rect_x_min <= x <= rect_x_max and rect_y_min <= y <= rect_y_max:
            overlap_count += 1

random.uniform() は呼び出されると 1 個の乱数しか返さない。そのため，必然的に返された乱数に基づく処理も 1 回ずつしか行えない，なので，for を使わざるを得なくなっている。

なお，処理部分はかなり複雑な計算をしているのに，結構な速さで for が回っているのだ。
num_points = 100000000 での実行速度は私の手元のマシンでは Wall time: 36.5 s であった。

例えば，これとは別の単純な（よくある正方形内の四分円で円周率を推計するような）やり方で for を使わずに numpy を使ってプログラムを書いてみると以下のようになる。

    x = np.random.rand(num_points)
    y = np.random.rand(num_points)
    overlap_count = np.sum((rect_x_min < x) * (x < rect_x_max) * (rect_y_min < y) * (y < rect_y_max) * ((x**2 + y**2) < 1))

このやり方は発生させた乱数対を少なからず無駄にするが，処理自体が単純なので前述のプログラムとほとんど同じ実行速度を持つ。
Wall time: 36.7 s であった。

逆に言えば，前述のプログラムを numpy を使って書くとずっと速くなるであろうということである。
以下のようになる。

    angle = np.random.rand(num_points)*2*math.pi
    r = np.sqrt(np.random.rand(num_points))*radius
    x = r * np.cos(angle)
    y = r * np.sin(angle)
    overlap_count = np.sum((rect_x_min <= x) * (x <= rect_x_max) * (rect_y_min <= y) * (y <= rect_y_max))

実行結果は Wall time: 4.12 s で，9 倍ほど速くなっている。

結論

numpy のベクトル計算をするには，道具自体も numpy.ndarray を対象にするような関数を使うことが必要だ。