arithmoi 0.1.0.0 → 0.1.0.1
raw patch · 4 files changed
+59/−26 lines, 4 filesPVP ok
version bump matches the API change (PVP)
API changes (from Hackage documentation)
Files
- Changes +4/−2
- Math/NumberTheory/Primes/Sieve.hs +34/−1
- Math/NumberTheory/Primes/Sieve/Eratosthenes.hs +20/−22
- arithmoi.cabal +1/−1
Changes view
@@ -1,2 +1,4 @@-0.1.0:-First release+0.1.0.1:+ Elaborate on overflow, work more on native Ints in Eratosthenes+0.1.0.0:+ First release
Math/NumberTheory/Primes/Sieve.hs view
@@ -16,7 +16,10 @@ -- where sieving is done, thus sieving primes up to @n@ requires -- @/O/(sqrt n/log n)@ space. module Math.NumberTheory.Primes.Sieve- ( primes+ ( -- * Limitations+ -- $limits+ -- * Sieves and lists+ primes , sieveFrom , PrimeSieve , primeSieve@@ -26,3 +29,33 @@ ) where import Math.NumberTheory.Primes.Sieve.Eratosthenes++-- $limits+--+-- There are three factors limiting the range of these sieves.+--+-- (1) Memory+--+-- (2) Overflow+--+-- (3) The internal representation of the state+--+-- An Eratosthenes type sieve needs to store the primes up to the square root of+-- the currently sieved region, thus requires @/O/(n\/log n)@ space.We store @16@ bytes+-- of information per prime, thus a Gigabyte of memory takes you to about @1.6*10^18@.+-- The @log@ doesn't change much in that range, so as a first approximation, doubling+-- the storage increases the sieve range by a factor of four.+--+-- On a 64-bit system, this is (currently) the only limitation to be concerned with, but+-- with more than four Terabyte of memory, the fact that the internal representation+-- currently limits the sieve range to about @6.8*10^25@ could become relevant.+-- Overflow in array indexing doesn't become a concern before memory and internal+-- representation would allow to sieve past @10^37@.+--+-- On a 32-bit system, the internal representation imposes no additional limits,+-- but overflow has to be reckoned with. On the one hand, the fact that arrays are+-- 'Int'-indexed restricts the size of the prime store, on the other hand, overflow+-- in calculating the indices to cross off multiples is possible before running out+-- of memory. The former limits the upper bound of the monolithic 'primeSieve' to+-- shortly above @8*10^9@, the latter limits the range of the segmented sieves to+-- about @1.7*10^18@.
Math/NumberTheory/Primes/Sieve/Eratosthenes.hs view
@@ -111,8 +111,9 @@ ] -- | List of primes.--- Since the sieve uses unboxed arrays, overflow occurs at some point,--- but not before @10^6*'fromIntegral' ('maxBound' :: 'Int')@ (I forgot where exactly).+-- Since the sieve uses unboxed arrays, overflow occurs at some point.+-- On 64-bit systems, that point is beyond the memory limits, on+-- 32-bit systems, it is at about @1.7*10^18@. primes :: [Integer] primes = 2:3:5:concat [[vO + toPrim i | i <- [0 .. li], unsafeAt bs i] | PS vO bs <- psieveList, let (_,li) = bounds bs]@@ -123,7 +124,7 @@ psieveList :: [PrimeSieve] psieveList = makeSieves plim sqlim 0 0 cache where- plim = 4801 -- prime #647+ plim = 4801 -- prime #647, 644 of them to use sqlim = plim*plim cache = runSTUArray $ do sieve <- sieveTo 4801@@ -135,12 +136,9 @@ if p then do let !i = indx .&. 7- k :: Integer- k = fromIntegral (indx `shiftR` 3)- strt1 = (k*(30*k + fromIntegral (2*rho i))- + fromIntegral (byte i)) `shiftL` 3- + fromIntegral (idx i)- !strt = fromIntegral strt1 .&. 0xFFFFF+ k = indx `shiftR` 3+ strt1 = (k*(30*k + 2*rho i) + byte i) `shiftL` 3 + fromIntegral (idx i)+ !strt = fromIntegral (strt1 .&. 0xFFFFF) !skip = fromIntegral (strt1 `shiftR` 20) !ixes = fromIntegral indx `shiftL` 23 + strt `shiftL` 3 + fromIntegral i unsafeWrite new j skip@@ -185,18 +183,18 @@ then unsafeWrite cache pr (w-1) else do ixes <- unsafeRead cache (pr+1)- let !stj = ixes .&. 0x7FFFFF- !ixw = ixes `shiftR` 23- !i = fromIntegral (ixw .&. 7)- !k = fromIntegral ixw - i+ let !stj = fromIntegral ixes .&. 0x7FFFFF -- position of multiple and index of cofactor+ !ixw = fromIntegral (ixes `shiftR` 23) -- prime data, up to 41 bits+ !i = ixw .&. 7+ !k = ixw - i -- On 32-bits, k > 44717396 means overflow is possible in tick !o = i `shiftL` 3- !j = fromIntegral stj .&. 7- !s = fromIntegral stj `shiftR` 3+ !j = stj .&. 7 -- index of cofactor+ !s = stj `shiftR` 3 -- index of first multiple to tick off (n, u) <- tick k o j s- let !skip = fromIntegral n `shiftR` 20- !strt = fromIntegral n .&. 0xFFFFF+ let !skip = fromIntegral (n `shiftR` 20)+ !strt = fromIntegral (n .&. 0xFFFFF) unsafeWrite cache pr skip- unsafeWrite cache (pr+1) (ixes - stj + strt `shiftL` 3 + fromIntegral u)+ unsafeWrite cache (pr+1) ((ixes .&. complement 0x7FFFFF) .|. strt `shiftL` 3 .|. fromIntegral u) treat (pr+2) tick stp off j ix | lastIndex < ix = return (ix - sieveBits, j)@@ -244,8 +242,8 @@ (bt,ix) = idxPr plim !start = 8*bt+ix+1 !nlim = plim+4800- sieve <- sieveTo nlim- (_,hi) <- getBounds sieve+ sieve <- sieveTo nlim -- Implement SieveFromTo for this, it's pretty wasteful when nlim isn't+ (_,hi) <- getBounds sieve -- very small anymore more <- countFromToWd start hi sieve new <- unsafeNewArray_ (0,num+2*more) :: ST s (STUArray s Int CacheWord) let copy i@@ -269,7 +267,7 @@ strt1 = strt0 - offset !strt = fromIntegral strt1 .&. 0xFFFFF !skip = fromIntegral (strt1 `shiftR` 20)- !ixes = fromIntegral indx `shiftL` 23 + strt `shiftL` 3 + fromIntegral i+ !ixes = fromIntegral indx `shiftL` 23 .|. strt `shiftL` 3 .|. fromIntegral i unsafeWrite new j skip unsafeWrite new (j+1) ixes fill (j+2) (indx+1)@@ -373,7 +371,7 @@ | otherwise = (strt1 - bitOff, i) !strt = fromIntegral strt2 .&. 0xFFFFF !skip = fromIntegral (strt2 `shiftR` 20)- !ixes = fromIntegral indx `shiftL` 23 + strt `shiftL` 3 + fromIntegral r2+ !ixes = fromIntegral indx `shiftL` 23 .|. strt `shiftL` 3 .|. fromIntegral r2 unsafeWrite new j skip unsafeWrite new (j+1) ixes fill (j+2) (indx+1)
arithmoi.cabal view
@@ -1,5 +1,5 @@ name : arithmoi-version : 0.1.0.0+version : 0.1.0.1 cabal-version : >= 1.6 author : Daniel Fischer copyright : (c) 2011 Daniel Fischer