[Haskell-cafe] A tale of Project Euler
Andrew Coppin
andrewcoppin at btinternet.com
Wed Nov 28 16:31:53 EST 2007
Sebastian Sylvan wrote:
> On Nov 28, 2007 12:12 PM, Kalman Noel <kalman.noel at bluebottle.com> wrote:
>
>> Sebastian Sylvan:
>>
>>> primes :: [Integer]
>>> primes = 2 : filter (null . primeFactors) [3,5..]
>>>
>>> primeFactors :: Integer-> [Integer]
>>> primeFactors n = factor n primes
>>> where
>>> factor m (p:ps) | p*p > m = []
>>> | m `mod` p == 0 = p : factor (m `div` p) (p:ps)
>>> | otherwise = factor m ps
>>>
>> Your definition gives a strange meaning to primeFactors. I'd want that for all
>> n, product (primeFactors n) == n. I think this law holds for the code posted
>> by Olivier.
>>
>
> Yes you're right. That is property should clearly hold.
>
There are problems for which it's important to know how many times a
given prime factor occurs. And there are other problems where it is
merely necessary to know which primes are factors. I would say it's
useful to have *both* functions available...
By the way, I'm now using a new and much uglier prime seive:
size = 1000000 :: Word64
calc_primes :: IO [Word64]
calc_primes = do
grid <- newArray (2,size) True
seive 2 grid
where
seive :: Word64 -> IOUArray Word64 Bool -> IO [Word64]
seive p g = do
mapM_ (\n -> writeArray g n False) [p, 2*p .. size]
mp' <- next (p+1) g
case mp' of
Nothing -> return [p]
Just p' -> do
ps <- seive p' g
return (p:ps)
next :: Word64 -> IOUArray Word64 Bool -> IO (Maybe Word64)
next p g = do
if p == size
then return Nothing
else do
t <- readArray g p
if t
then return (Just p)
else next (p+1) g
Seems moderately fast. At least, I can find the 10,001st prime in under
1 second... The obvious downside of course is that you must know how big
to make the array at the start of your program, and you must compute the
entire array before you can use any of it. Both limitations not precent
in the lazy recursive list implementation...
More information about the Haskell-Cafe
mailing list