MarcoXerox
MarcoXerox

Reputation: 68

Haskell vs C speed: Sieve of Eratosthenes

Haskell one is implemented using optimized Data.IntSet with complexity O(lg n). However, there is a 15x (previously 30x) speed difference for n = 2000000 despite Haskell code is already optimized for even number cases. I would like to know whether/why my implementation in Haskell is imperfect.

Original Haskell

primesUpTo :: Int -> [Int]
primesUpTo n = 2 : put S.empty [3,5..n]
    where put :: S.IntSet -> [Int] -> [Int]
          put _ [] = []
          put comps (x:xs) =
            if S.member x comps
            then put comps xs
            else x : put (S.union comps multiples) xs
                where multiples = S.fromList [x*2, x*3 .. n]

Update

fromDistinctAscList gives a 4x speed increase. 2-3-5-7-Wheel speeds up by another 50%.

primesUpTo :: Int -> [Int]
primesUpTo n = 2 : 3 : 5 : 7 : put S.empty (takeWhile (<=n) (spin wheel 11))
    where put :: S.IntSet -> [Int] -> [Int]
          put _ [] = []
          put comps (x:xs) =
            if S.member x comps
            then put comps xs
            else x : put (S.union comps multiples) xs
                where multiples = S.fromDistinctAscList [x*x, x*(x+2) .. n]
          spin (x:xs) n = n : spin xs (n + x)
          wheel = 2:4:2:4:6:2:6:4:2:4:6:6:2:6:4:2:6:4:6:8:4:2:4:2:4:8:6:4:6:2:4:6:2:6:6:4:2:4:6:2:6:4:2:4:2:10:2:10:wheel

Benchmarking

All time are measured by *nix time command, real space

Haskell original : 2e6: N/A;    2e7: >30s
Haskell optimized: 2e6: 0.396s; 2e7: 6.273s
C++ Set (ordered): 2e6: 4.694s; 2e7: >30s
C++ Bool Array   : 2e6: 0.039s; 2e7: 0.421s

Haskell optimized is slower than C++ Bool by 10~15x, and faster than C++ Set by 10x.

Source code

C Compiler options: g++ 5.3.1, g++ -std=c++11 Haskell options: ghc 7.8.4, ghc

C code (Bool array) http://pastebin.com/W0s7cSWi

 prime[0] = prime[1] = false;
 for (int i=2; i<=limit; i++) { //edited
     if (!prime[i]) continue;
     for (int j=2*i; j<=n; j+=i)
        prime[j] = false;
 }

C code (Set) http://pastebin.com/sNpghrU4

 nonprime.insert(1);
 for (int i=2; i<=limit; i++) { //edited
     if (nonprime.count(i) > 0) continue;
     for (int j=2*i; j<=n; j+=i)
        nonprime.insert(j);
 }

Haskell code http://pastebin.com/HuMqwvRW Code as written above.

Upvotes: 2

Views: 1115

Answers (1)

behzad.nouri
behzad.nouri

Reputation: 78021

I would like to know whether/why my implementation in Haskell is imperfect.

Instead of fromList you better use fromDistinctAscList which performs linearly. You may also add only odd multiples starting with x*x not x*2, because all the smaller odd multiples have already been added. Style-wise, a right fold may fit better than recursion.

Doing so, I get more than 3 times performance improvement for n equal to 2,000,000:

import Data.IntSet (member, union, empty, fromDistinctAscList)

sieve :: Int -> [Int]
sieve n = 2: foldr go (const []) [3,5..n] empty
    where
    go i run obs
        | member i obs = run obs
        | otherwise    = i: run (union obs inc)
        where inc = fromDistinctAscList [i*i, i*(i + 2)..n]

Nevertheless, an array has both O(1) access and cache friendly memory allocation. Using mutable arrays, I see more than 15 times performance improvement over your Haskell code (again n equal to 2,000,000):

{-# LANGUAGE FlexibleContexts #-}
import Data.Array.ST (STUArray)
import Control.Monad (forM_, foldM)
import Control.Monad.ST (ST, runST)
import Data.Array.Base (newArray, unsafeWrite, unsafeRead)

sieve :: Int -> [Int]
sieve n = reverse $ runST $ do
    arr <- newArray (0, n) False :: ST s (STUArray s Int Bool)
    foldM (go arr) [2] [3,5..n]
    where
    go arr acc i = do
        b <- unsafeRead arr i
        if b then return acc else do
            forM_ [i*i, i*(i + 2).. n] $ \k -> unsafeWrite arr k True
            return $ i: acc

Upvotes: 4

Related Questions