[Go] 並行処理を用いたモンテカルロ法

Go言語による並行処理の第4章4.6パイプラインをもとにモンテカルロ法で円周率を計算する。 www.amazon.co.jp

乱数による円周率の計算

以下の積分モンテカルロ法を適用し、円周率をモンテカルロ法で求める。 $$ \int_{0}^{1} \frac{4}{ 1 + x^{2} } dx = \pi $$

サンプル数 $ M $ を$ 10,10^{3},10^{5},10^{7} $として、

$$ \bar{f_i} = \frac{1}{M} \sum_{t = 1, ..., M} \frac{4}{1 + {X_t}^{2}} $$

を求め、推定値

$$ \bar{ \bar{f} } = \frac{\sum_i \bar{f}_i}{N_M} $$

を計算する。ただし、$ N_{M}$は独立なモンテカルロ計算の回数である。

並行処理

パイプライン: データを受け取って、何らかの処理をし、どこかに渡すという一連の処理。

本に書いてあるrepeatFn関数やtake関数を用いて実装する。

実装

func repeatFn(done <-chan interface{}, fn func() float64) <-chan float64 {
    valueStream := make(chan float64)
    go func() {
        defer close(valueStream)
        for {
            select {
            case <-done:
                return
            case valueStream <- fn():
            }
        }
    }()
    return valueStream
}

func take(done <-chan interface{}, valueStream <-chan float64, num int) <-chan float64 {
    takeStream := make(chan float64)
    go func() {
        defer close(takeStream)
        for i := 0; i < num; i++ {
            select {
            case <-done:
                return
            case takeStream <- <-valueStream:
            }
        }
    }()
    return takeStream
}

func calc(done <-chan interface{}, valueStream <-chan float64, fn func(float64) float64) <-chan float64 {
    calcStream := make(chan float64)
    go func() {
        defer close(calcStream)
        for x := range valueStream {
            select {
            case <-done:
                return
            case calcStream <- fn(x):
            }
        }
    }()
    return calcStream
}

func calculatePi(n int) float64 {
    var piApprox float64
    for i := 0; i < n; i++ {
        x := rand.Float64()
        piApprox += (4 / (1 + x*x))
    }
    piApprox /= float64(n)
    return piApprox
}

func main() {
    done := make(chan interface{})
    defer close(done)

    var piApprox float64
    var n float64 = math.Pow(10, 6)
    rand := func() float64 { return rand.Float64() }
    f := func(x float64) float64 { return 4 / ((1 + x*x) * n) }

    start := time.Now()
    loop := func() <-chan float64 {
        valueStream := make(chan float64)
        go func() {
            defer close(valueStream)
            for i := 0; i < 10; i++ {
                var mean float64
                for num := range calc(done, take(done, repeatFn(done, rand), int(n)), f) {
                    mean += num
                }
                valueStream <- mean
            }
        }()
        return valueStream
    }

    for num := range loop() {
        piApprox += num
    }

    fmt.Printf("calculated pi: %v, time: %v s\n", piApprox/10, time.Since(start).Seconds())

    start = time.Now()

    piApprox = 0
    for i := 0; i < 10; i++ {
        piApprox += calculatePi(int(n))
    }
    fmt.Printf("calculated pi: %v, time: %v s\n", piApprox/10, time.Since(start).Seconds())

}

calc関数では乱数$ x$に対して$ \frac{4}{1 + x^{2}}$を計算し、それをcalcStreamに追加している。

また、calculatePiは並行処理をせずに単純にモンテカルロ計算をするための関数である。

結果

  • 並行処理
サンプル数 計算された円周率 実行時間(s)
$10$ $ 3.15711$ $ 4.36634 \times 10^{-4}$
$10^{3}$ $ 3.14847$ $0.018$
$10^{5}$ $ 3.14215$ $1.771$
$10^{7}$ $ 3.14151$ $174.374$
  • 単純に処理
サンプル数 計算された円周率 実行時間(s)
$10$ $ 3.23199$ $2.757 \times 10^{-6}$
$10^{3}$ $ 3.14217$ $ 1.75575 \times 10^{-4}$
$10^{5}$ $ 3.14065$ $0.019$
$10^{7}$ $ 3.141571$ $1.772$

考察

明らかに単純にcalculatePi関数を用いて計算した方が速いです。何重にもするのではなく、もっとパイプラインをまとめたほうが速くなりそうです。

また、パイプラインの使い方を本を元に実装する上でちょうど良い例になったのではないかと思います。