前述のように、可変配列を使用すると最高のパフォーマンスが得られます。次のコードは、この「TemplateHaskell」バージョンから派生したもので、「TemplateHaskell」は何の違いもないように見えるため、さらにいくつかの最適化を加えて、単純な可変配列ソリューションに沿ったものに変換し直しています。さらなる最適化と、特に配列の範囲チェックを回避する「unsafeRead」および「unsafeWrite」関数の使用により、通常の可変ボックス化されていない配列バージョンよりも高速であると思います。おそらく配列アクセス用のポインターも内部的に使用します。
{-# OPTIONS -O2 -optc-O3 #-}
import Control.Monad
import Control.Monad.ST
import Data.Array.ST
import Data.Array.Unboxed
import Data.Array.Base
primesToUA :: Word32-> [Word32]
primesToUA top = do
let sieveUA top = runSTUArray $ do
let m = ((fromIntegral top) - 3) `div` 2 :: Int
buf <- newArray (0,m) True -- :: ST s (STUArray s Int Bool)
let cullUA i = do
let p = i + i + 3
strt = p * (i + 1) + i
let cull j = do
if j > m then cullUA (i + 1)
else do
unsafeWrite buf j False
cull (j + p)
if strt > m then return ()
else do
e <- unsafeRead buf i
if e then cull strt else cullUA (i + 1)
cullUA 0; return buf
if top > 1 then 2 : [2 * (fromIntegral i) + 3 | (i,True) <- assocs $ sieveUA top]
else []
main = do
x <- read `fmap` getLine -- 1mln 2mln 10mln 100mln
print (last (primesToUA x)) -- 0.01 0.02 0.09 1.26 seconds
編集: 上記のコードが修正され、以下の時間が修正を反映するように編集されました。また、ページ付きバージョンとページなしバージョンを比較するコメントも編集されました。
これを指定された上限範囲まで実行する時間は、上記のソース コードの下部にあるコメント テーブルに示されているとおりであり、ideone.com によって測定され、 @WillNess によって投稿された回答 ( ideoneで測定されたもの) よりも正確に約 5 倍高速です。 .com. このコードは、素数を 200 万に選別するのにわずかな時間を要し、1 億に選別するのにわずか 1.26 秒しかかかりません。私の i7 (3.5 GHz) で 0.44 秒から 1 億まで実行すると、これらの時間は約 2.86 倍速くなり、10 億まで実行するには 6.81 秒かかります。メモリ使用量は、前者で 6 メガバイト強、後者で 60 メガバイト強です。これは、巨大な (ビット パックされた) 配列によって使用されるメモリです。この配列はまた、配列サイズが CPU キャッシュ サイズを超えると、合成数表現カリングごとに平均メモリ アクセス時間が悪化するという点で、非線形パフォーマンスを説明します。
編集_追加: ページ分割されたふるいは、バッファー サイズが L1 または L2 CPU キャッシュよりも小さく保たれている場合にメモリ アクセス効率が向上するという点でより効率的であり、上限がなくてもよいように無制限であるという利点もあります。事前に指定され、使用される範囲の平方根にページ バッファ メモリを加えた値よりも小さい基本素数だけである、はるかに小さいメモリ フットプリントです。次のコードは、ページ セグメント化された実装として記述されており、ページ化されていないバージョンよりも多少高速です。また、コードの先頭にある出力範囲の指定を「Word32」から「Word64」に変更できるという利点もあります。そのため、処理時間のわずかなコストで、32 ビットの数値範囲に限定されません ( 32 ビットのコンパイル済みコードの場合) 共通の任意の範囲。コードは次のとおりです。
-- from http://www.haskell.org/haskellwiki/Prime_numbers#Using_ST_Array
{-# OPTIONS -O2 -optc-O3 #-}
import Data.Word
import Control.Monad
import Control.Monad.ST
import Data.Array.ST
import Data.Array.Unboxed
import Data.Array.Base
primesUA :: () -> [Word32]
primesUA () = do
let pgSZBTS = 262144 * 8
let sieveUA low bps = runSTUArray $ do
let nxt = (fromIntegral low) + (fromIntegral pgSZBTS)
buf <- newArray (0,pgSZBTS - 1) True -- :: ST s (STUArray s Int Bool)
let cullUAbase i = do
let p = i + i + 3
strt = p * (i + 1) + i
when (strt < pgSZBTS) $ do
e <- unsafeRead buf i
if e then do
let loop j = do
if j < pgSZBTS then do
unsafeWrite buf j False
loop (j + p)
else cullUAbase (i + 1)
loop strt
else cullUAbase (i + 1)
let cullUA ~(p:t) = do
let bp = (fromIntegral p)
i = (bp - 3) `div` 2
s = bp * (i + 1) + i
when (s < nxt) $ do
let strt = do
if s >= low then fromIntegral (s - low)
else do
let b = (low - s) `rem` bp
if b == 0 then 0 else fromIntegral (bp - b)
let loop j = do
if j < pgSZBTS then do
unsafeWrite buf j False
loop (j + (fromIntegral p))
else cullUA t
loop strt
if low <= 0 then cullUAbase 0 else cullUA bps
return buf
let sieveList low bps = do
[2 * ((fromIntegral i) + low) + 3 | (i,True) <- assocs $ sieveUA low bps]
let sieve low bps = do
(sieveList low bps) ++ sieve (low + (fromIntegral pgSZBTS)) bps
let primes' = ((sieveList 0 []) ++ sieve (fromIntegral pgSZBTS) primes') :: [Word32]
2 : sieve 0 primes'
main = do
x <- read `fmap` getLine -- 1mln 2mln 10mln 100mln
-- 0.02 0.03 0.13 1.13 seconds
print (length (takeWhile ((>=) (fromIntegral x)) (primesUA ())))
上記のコードは、後続のページとは異なる方法で最初のページ配列から合成数表現を選別する必要があるため、非ページの場合よりもかなり多くのコード行があります。コードには修正も含まれているため、基本素数リストと出力リストが同じリストではないためにメモリ リークが発生することはありません (したがって、メモリ内のリスト全体を保持することを回避できます)。
このコードは、カリング バッファーが L2 CPU キャッシュよりも小さい一定サイズであるために範囲が大きくなるため、(ふるいにかけられた範囲にわたって) 線形時間に近くなることに注意してください。メモリ使用量は、1 億の場合は 600 キロバイト弱、10 億の場合は 600 キロバイト強で、非ページ バージョンで使用されるものの一部です。わずかな増加は、平方根未満の基本素数に必要な余分なスペースにすぎません。範囲リストの。
ideone.comでは、このコードにより、約 1.13 秒で 1 億個までの素数が生成され、10 億個までは約 12 秒で生成されます (32 ビット設定)。おそらくホイールの因数分解と、間違いなくマルチコア処理により、マルチコア CPU でさらに高速になります。私の i7 (3.5 GHz) では、1 億にふるいにかけるのに 0.44 秒、10 億にふるいにかけるのに 4.7 秒かかります。ideone.com で実行されている GHC のバージョンには、i7 には存在しない、より大きな範囲でパフォーマンスが低下する、ある種の非線形オーバーヘッドがあるようです。新しいページごとにバッファが新しく作成されています。 END_EDIT_ADD
EDIT_ADD2:上記のページ セグメント化されたコードの処理時間の多くは (遅延) リスト処理で使用されているようです。したがって、コードは次のようにいくつかの改善を加えて再編成されています。
リスト処理を使用せず、「popCount」テーブル ルックアップを使用して一度に 32 ビット ワード内の「1」ビットの数をカウントする素数カウント関数を実装しました。このように、結果を見つけるための時間は、実際のふるい選別時間と比較して重要ではありません。
基本素数をビット パックされたページ セグメントのリストとして格納しました。これは、素数のリストを格納するよりもはるかにスペース効率が高く、必要に応じてページ セグメントを素数に変換する時間は計算オーバーヘッドが大きくありません。
最初のゼロ ページ セグメントに対して、独自のビット パターンをソース ページとして使用するように素数セグメントの make 関数を調整し、合成数カリング コードをより短く、より単純にしました。
コードは次のようになります。
{-# OPTIONS -O3 -rtsopts #-} -- -fllvm ide.com doesn't support LLVM
import Data.Word
import Data.Bits
import Control.Monad
import Control.Monad.ST
import Data.Array.ST (runSTUArray)
import Data.Array.Unboxed
import Data.Array.Base
pgSZBTS = (2^18) * 8 :: Int -- size of L2 data cache
type PrimeType = Word32
type Chunk = UArray PrimeType Bool
-- makes a new page chunk and culls it
-- if the base primes list provided is empty then
-- it uses the current array as source (for zero page base primes)
mkChnk :: Word32 -> [Chunk] -> Chunk
mkChnk low bschnks = runSTUArray $ do
let nxt = (fromIntegral low) + (fromIntegral pgSZBTS)
buf <- nxt `seq` newArray (fromIntegral low, fromIntegral nxt - 1) True
let cull ~(p:ps) =
let bp = (fromIntegral p)
i = (bp - 3) `shiftR` 1
s = bp * (i + 1) + i in
let cullp j = do
if j >= pgSZBTS then cull ps
else do
unsafeWrite buf j False
cullp (j + (fromIntegral p)) in
when (s < nxt) $ do
let strt = do
if s >= low then fromIntegral (s - low)
else do
let b = (low - s) `rem` bp
if b == 0 then 0 else fromIntegral (bp - b)
cullp strt
case bschnks of
[] -> do bsbf <- unsafeFreezeSTUArray buf
cull (listChnkPrms [bsbf])
_ -> cull $ listChnkPrms bschnks
return buf
-- creates a page chunk list starting at the lw value
chnksList :: Word32 -> [Chunk]
chnksList lw =
mkChnk lw basePrmChnks : chnksList (lw + fromIntegral pgSZBTS)
-- converts a page chunk list to a list of primes
listChnkPrms :: [Chunk] -> [PrimeType]
listChnkPrms [] = []
listChnkPrms ~(hdchnk@(UArray lw _ rng _):tlchnks) =
let nxtp i =
if i >= rng then [] else
if unsafeAt hdchnk i then
(case ((lw + fromIntegral i) `shiftL` 1) + 3 of
np -> np) : nxtp (i + 1)
else nxtp (i + 1) in
(hdchnk `seq` lw `seq` nxtp 0) ++ listChnkPrms tlchnks
-- the base page chunk list used to cull the higher order pages,
-- note that it has special treatment for the zero page.
-- It is more space efficient to store this as chunks rather than
-- as a list of primes or even a list of deltas (gaps), with the
-- extra processing to convert as needed not too much.
basePrmChnks :: [Chunk]
basePrmChnks = mkChnk 0 [] : chnksList (fromIntegral pgSZBTS)
-- the full list of primes could be accessed with the following function.
primes :: () -> [PrimeType]
primes () = 2 : (listChnkPrms $ chnksList 0)
-- a quite fast prime counting up to the given limit using
-- chunk processing to avoid innermost list processing loops.
countPrimesTo :: PrimeType -> Int
countPrimesTo limit =
let lmtb = (limit - 3) `div` 2 in
let sumChnks acc chnks@(chnk@(UArray lo hi rng _):chnks') =
let cnt :: UArray PrimeType Word32 -> Int
cnt bfw =
case if lmtb < hi then fromIntegral (lmtb - lo) else rng of
crng -> case crng `shiftR` 5 of
rngw ->
let cnt' i ac =
ac `seq` if i >= rngw then
if (i `shiftL` 5) >= rng then ac else
case (-2) `shiftL` fromIntegral (lmtb .&. 31) of
msk -> msk `seq`
case (unsafeAt bfw rngw) .&.
(complement msk) of
bts -> bts `seq` case popCount bts of
c -> c `seq` case ac + c of nac -> nac
else case ac + (popCount $ unsafeAt bfw i) of
nacc -> nacc `seq` cnt' (i + 1) (nacc)
in cnt' 0 0 in
acc `seq` case runST $ do -- make UArray _ Bool into a UArray _ Word32
stbuf <- unsafeThawSTUArray chnk
stbufw <- castSTUArray stbuf
bufw <- unsafeFreezeSTUArray stbufw
return $ cnt bufw of
c -> c `seq` case acc + c of
nacc -> nacc `seq` if hi >= lmtb then nacc
else sumChnks nacc chnks' in
if limit < 2 then 0 else if limit < 3 then 1 else
lmtb `seq` sumChnks 1 (chnksList 0)
main = do
x <- read `fmap` getLine -- 1mln 2mln 10mln 100mln 1000mln
-- 0.02 0.03 0.06 0.45 4.60 seconds
-- 7328 7328 8352 8352 9424 Kilobytes
-- this takes 14.34 seconds and 9424 Kilobytes to 3 billion on ideone.com,
-- and 9.12 seconds for 3 billion on an i7-2700K (3.5 GHz).
-- The above ratio of about 1.6 is much better than usual due to
-- the extremely low memory use of the page segmented algorithm.
-- It seems thaat the Windows Native Code Generator (NCG) from GHC
-- is particularly poor, as the Linux 32-bit version takes
-- less than two thirds of the time for exactly the same program...
print $ countPrimesTo x
-- print $ length $ takeWhile ((>=) x) $ primes () -- the slow way to do this
コードに示されている時間とメモリ要件は、ideone.com で実行したときに観察されたものです。、0.02、0.03、0.05、0.30、3.0、および 9.1 秒で、i7-2700K (3.5 GHz) で 1、2、10、100、1000 (10 億)、および 3000 (30 億) を実行する必要があります。 ) 100 万の範囲であり、メモリ フットプリントはほぼ一定であり、必要に応じて範囲の平方根よりも少ない基本素数の数と共にゆっくりと増加します。LLVM コンパイラ バックエンドでコンパイルすると、これらの時間はそれぞれ 0.01、0.02、0.02、0.12、1.35、および 4.15 秒になります。これは、レジスタとマシン命令をより効率的に使用するためです。この最後のものは、レジスターを効率的に使用するために使用される 32 ビット コンパイラーではなく、64 ビット コンパイラーでコンパイルされた場合と同じ速度に非常に近く、余分なレジスターの可用性は大きな違いをもたらさないことを意味します。
コードでコメントされているように、私の実機と ideone.com サーバーのパフォーマンスの比率は、メモリ アクセスのボトルネックによって抑制されないため、メモリを無駄に消費するアルゴリズムよりもはるかに小さくなり、速度の制限はほとんどの場合、 CPU クロック速度とクロック サイクルあたりの CPU 処理効率。ただし、そこでコメントされているように、Windows (32 ビット コンパイラ) で GHC ネイティブ コード ジェネレーター (NCG) を実行すると、奇妙な非効率性があり、Linux で実行する場合よりも実行時間が 50% 以上遅くなります (イデオンとして)。 com サーバーが使用します)。私の知る限り、どちらも同じHaskell GHCバージョンの共通コードベースを持っており、GHC NCGはGCCを使用せず、mingw32アセンブラーのみを使用するため、唯一の違いは使用されるリンカーにあります(これは影響を受けないLLVMバックエンドでも使用されます)。 、
このコードを LLVM コンパイラ バックエンドでコンパイルすると、高度に最適化された「C/C++」実装用に記述された同じアルゴリズムとほぼ同じ速度になることに注意してください。これは、Haskell が非常にタイトなループ コーディングを開発する能力を実際に備えていることを示しています。Haskell のコードは、Monadic および Non-strict コードの Haskell パラダイムに慣れると、同等の「C/C++」コードよりもかなり読みやすく、安全であると言えます。エラトステネスの篩の実行速度のさらなる改善は、純粋に使用される実装の調整の関数であり、Haskell と 'C/C++' の間の言語の選択ではありません。
要約: もちろん、これはエラトステネスのふるいの Haskell バージョンとしてはまだ究極の速度ではありません。つまり、高速な CPU L1 キャッシュをより効率的に使用するためのメモリ アクセスをまだ調整しておらず、パフォーマンスを大幅に低下させていません。オッズ処理を排除する以外に極端なホイール因数分解を使用して必要な複合カリング操作の総数。ただし、リストや不変配列を使用するよりも約 100 倍高速になる可能性があり、可変配列がこのようなタイト ループ タイプの問題に対処する最も効率的な方法であることを示すには、これで十分です。 END_EDIT_ADD2