• 作成:

Re:Haskellで書いてみたらC++の10倍遅かった 5倍程度になりました

10倍は遅すぎませんか?

本の虫: Haskellで書いてみたらC++の10倍遅かった

を読みました.

10倍は差が出過ぎなのではないかと思いました.

C++ソースコード

まずC++のソースコードが完全なものではないので補完しました.

#include <algorithm>
#include <array>
#include <iostream>
#include <random>

template <typename Random>
bool coin_toss(const unsigned int try_count, const unsigned int length,
               Random &r) {
  unsigned int count{};
  int prev{};

  std::uniform_int_distribution<> d(0, 1);

  for (unsigned int i = 0; i != try_count; ++i) {
    int result = d(r);
    if (result == prev) {
      ++count;
      if (count == length)
        return true;
    } else {
      prev = result;
      count = 1;
    }
  }

  return false;
}

template <typename Random>
double monte_carlo(unsigned int try_count, Random &r) {
  unsigned int count{};
  for (unsigned int i = 0; i != try_count; ++i) {
    if (coin_toss(100, 10, r))
      ++count;
  }

  return double(count) / double(try_count);
}

int main() {
  std::array<unsigned int, sizeof(std::mt19937) / sizeof(unsigned int)> c;
  std::generate(std::begin(c), std::end(c), std::random_device{});
  std::seed_seq s(std::begin(c), std::end(c));
  std::mt19937 engine(s);

  for (unsigned int i = 100; i != 1000000; i *= 10) {
    auto r = engine;
    std::cout << i << "\t : " << monte_carlo(i, r) * 100.0 << "%\n";
  }
}
2018-02-06T12:05:17 ncaq@strawberry/pts/4(1) ~/Documents/archive/2018-02
% g++ -O2 -std=c++17 coin.cpp
2018-02-06T12:05:30 ncaq@strawberry/pts/4(0) ~/Documents/archive/2018-02
% time ./a.out
100	 : 14%
1000	 : 8.3%
10000	 : 8.8%
100000	 : 8.557%
./a.out  0.21s user 0.00s system 99% cpu 0.217 total

Haskellソースコードにhlintをかけたもの

Haskellソースコードを写して, 元のソースコードが読みにくいのでhlintに従って整形しました.

import           Data.List
import           System.Random

coinSeq :: (RandomGen g) => g -> [Bool]
coinSeq = randoms

splitN :: Int -> [a] -> [[a]]
splitN n s = take n s : splitN n (drop n s)

hasContiguousElems :: (Eq a) => Int -> [a] -> Bool
hasContiguousElems len s = any (>= len) (map length (group s))

coinToss :: Int -> [Bool] -> Bool
coinToss = hasContiguousElems

monteCarlo :: Int -> [Bool] -> Double
monteCarlo try_count s = (fromIntegral n / fromIntegral try_count) * 100.0
    where
        seq_n = take try_count (splitN 100 s)
        n = length . filter id $ map (coinToss 10) seq_n

main :: IO ()
main = do
    gen <- getStdGen
    let s = coinSeq gen
    mapM_ (\n -> putStrLn (show n ++ "\t : " ++ show (monteCarlo n s) ++ "%") )
        [100, 1000, 10000, 100000]
2018-02-06T13:43:57 ncaq@strawberry/pts/4(130) ~/Documents/archive/2018-02
% stack ghc -- -O2 oldCoin.hs
2018-02-06T13:43:59 ncaq@strawberry/pts/4(0) ~/Documents/archive/2018-02
% time ./oldCoin
100	 : 10.0%
1000	 : 8.1%
10000	 : 8.780000000000001%
100000	 : 8.559999999999999%
./oldCoin  3.23s user 0.03s system 99% cpu 3.268 total

なるほど確かに10倍程度遅いようですね.

プロファイリング

ところでghcにはプロファイリングオプションがあります.

計測せよって偉い人も言ってました.

2018-02-06T12:11:55 ncaq@strawberry/pts/6(0) ~/Documents/archive/2018-02
% stack ghc -- -prof -fprof-auto -rtsopts coin.hs
[1 of 1] Compiling Main             ( coin.hs, coin.o )
Linking coin ...
2018-02-06T12:13:11 ncaq@strawberry/pts/6(0) ~/Documents/archive/2018-02
% time ./coin +RTS -p
100	 : 91.0%
1000	 : 91.2%
10000	 : 91.43%
100000	 : 91.259%
./coin +RTS -p  10.76s user 0.11s system 99% cpu 10.882 total
        Tue Feb  6 12:13 2018 Time and Allocation Profiling Report  (Final)

           coin +RTS -p -RTS

        total time  =       10.15 secs   (10151 ticks @ 1000 us, 1 processor)
        total alloc = 11,164,362,552 bytes  (excludes profiling overheads)

COST CENTRE               MODULE        SRC                                %time %alloc

hasContiguousElems        Main          coin.hs:11:1-68                     21.2   25.6
randomIvalInteger         System.Random System/Random.hs:(468,1)-(489,76)   21.0   26.5
stdNext                   System.Random System/Random.hs:(518,1)-(528,64)   13.8   15.8
randomIvalInteger.f       System.Random System/Random.hs:(486,8)-(489,76)    9.3    2.1
randomIvalInteger.f.v'    System.Random System/Random.hs:489:25-76           5.3    4.3
randomIvalInteger.b       System.Random System/Random.hs:473:8-54            4.0    5.7
randomIvalInteger.(...)   System.Random System/Random.hs:472:8-36            3.8    0.0
randomIvalInteger.magtgt  System.Random System/Random.hs:483:8-21            3.7    1.4
randomR                   System.Random System/Random.hs:(387,3)-(397,30)    2.9    0.0
next                      System.Random System/Random.hs:218:3-17            2.8    2.1
splitN                    Main          coin.hs:8:1-43                       2.7    5.4
randoms.\                 System.Random System/Random.hs:316:42-67           2.4    5.0
randomIvalInteger.f.(...) System.Random System/Random.hs:488:25-39           2.3    0.0
randomIvalInteger.k       System.Random System/Random.hs:482:8-20            1.8    1.4
random                    System.Random System/Random.hs:399:3-49            1.1    0.0
stdNext.z                 System.Random System/Random.hs:520:17-34           0.6    1.4
stdNext.s1''              System.Random System/Random.hs:524:17-64           0.3    1.4
stdNext.s2''              System.Random System/Random.hs:528:17-64           0.2    1.4

やっぱり乱数生成が遅い

あーやっぱりRandomが遅い. そんな気はしていました.

Haskellの乱数事情 - Qiita で触れられてるとおり, System.Randomは遅いのです.

乱数を速度を重視して真面目に考えたことがないので何を選択すれば良いのか全くわからない. 先程の記事で推奨されていたmwc-randomを使ってみます.

一番速いものを選ぼうかと思ったのですが, 純Haskellで実装されたものでないと「Haskellの速度」を測るレギュレーションに違反しているのではないかと思ったので, 純Haskellで実装されたものを選びました.

System.Random.MWCでの再実装

import           Data.List
import qualified Data.Vector       as V
import           System.Random.MWC

coinSeq :: Int -> IO (V.Vector Bool)
coinSeq l = withSystemRandom . asGenST $ \gen -> uniformVector gen l

splitN :: Int -> [a] -> [[a]]
splitN n s = take n s : splitN n (drop n s)

hasContiguousElems :: (Eq a) => Int -> [a] -> Bool
hasContiguousElems len s = any (>= len) (map length (group s))

coinToss :: Int -> [Bool] -> Bool
coinToss = hasContiguousElems

monteCarlo :: Int -> [Bool] -> Double
monteCarlo try_count s = (fromIntegral n / fromIntegral try_count) * 100.0
    where
        seq_n = take try_count (splitN 100 s)
        n = length . filter id $ map (coinToss 10) seq_n

main :: IO ()
main = do
    s <- V.toList <$> coinSeq 100000
    mapM_ (\n -> putStrLn (show n ++ "\t : " ++ show (monteCarlo n s) ++ "%") )
        [100, 1000, 10000, 100000]
2018-02-06T13:08:09 ncaq@strawberry/pts/6(0) ~/Documents/archive/2018-02
% stack ghc -- -O2 coin.hs
2018-02-06T13:08:13 ncaq@strawberry/pts/6(0) ~/Documents/archive/2018-02
% time ./coin
100	 : 9.0%
1000	 : 9.2%
10000	 : 0.9199999999999999%
100000	 : 9.2e-2%
./coin  0.03s user 0.00s system 98% cpu 0.038 total

これでC++の2倍程度の遅さに収まりましたね. 良かったですね.

と言いたいところですが, 乱数が偏りすぎていてまともな答えが帰ってきていない. これではコードがちゃんと動いたとは言い難いです.

まともな偏りを持った純Haskell製乱数生成ライブラリは少しだけ調べた結果見つかりませんでした. 実用的には他言語へのバインディングライブラリ使えば良いんでしょうけど, 今回のケースを解決するものは見つかりませんでした. 残念です.

答えとしては片手落ちですね. 誰か知ってる人が居たら教えてください.

そもそも無限リストなんかじゃなくてRandomIOとか使って必要なところに乱数を使えばそこそこ早くなるんじゃないのと思ったのですがとても眠いので大幅な書き換えをやる気が起きません.

指摘してもらいました

なるほど.

元のアルゴリズムを全く読まずに適当な書き換えをしたツケですね.

乱数に偏りがあったという私の理解は誤りでした.

修正版

import           Data.List
import qualified Data.Vector       as V
import           System.Random.MWC

coinSeq :: Int -> IO (V.Vector Bool)
coinSeq l = withSystemRandom . asGenST $ \gen -> uniformVector gen l

splitN :: Int -> [a] -> [[a]]
splitN n s = take n s : splitN n (drop n s)

hasContiguousElems :: (Eq a) => Int -> [a] -> Bool
hasContiguousElems len s = any (>= len) (map length (group s))

coinToss :: Int -> [Bool] -> Bool
coinToss = hasContiguousElems

monteCarlo :: Int -> [Bool] -> Double
monteCarlo try_count s = (fromIntegral n / fromIntegral try_count) * 100.0
    where
        seq_n = take try_count (splitN 100 s)
        n = length . filter id $ map (coinToss 10) seq_n

main :: IO ()
main = do
    s <- V.toList <$> coinSeq (100000 * 100)
    mapM_ (\n -> putStrLn (show n ++ "\t : " ++ show (monteCarlo n s) ++ "%") )
        [100, 1000, 10000, 100000]
2018-02-06T16:54:40 ncaq@strawberry/pts/9(0) ~/Documents/archive/2018-02
% time ./coin
100	 : 8.0%
1000	 : 7.7%
10000	 : 8.959999999999999%
100000	 : 8.698%
./coin  1.09s user 0.04s system 99% cpu 1.133 total

これで1.133 / 0.217 = 5.221198156682028というわけで, 5倍程度の速さ(遅さ)になりました.

リストじゃなくてVectorのまま取り扱ったらもっと早くなるだろうって? それは読者の宿題としようと思います(めんどくさい).