Posts Tagged ‘euler’

Programming Praxis – Sieve Of Euler

February 25, 2011

In today’s Programming Praxis exercise, our goal is to implement the Sieve of Euler to generate primes. Let’s get started, shall we?

Some imports:

import Data.List.Ordered
import Data.Time.Clock

Euler’s sieve is pretty simple: just filter out all remaining multiples of each successive prime number. We add the minor optimization of only considering odd numbers.

primes_euler :: [Integer]
primes_euler = 2 : euler [3,5..] where
    euler ~(x:xs) = x : euler (minus xs $ map (* x) (x:xs))

To compare, we use the Sieve of Eratosthenes implementation found here:

primes_eratosthenes :: [Integer]
primes_eratosthenes = 2 : 3 : sieve (tail primes_eratosthenes) 3 [] where
    sieve ~(p:ps) x fs = let q=p*p in
        foldr (flip minus) [x+2,x+4..q-2] [[y+s,y+2*s..q] | (s,y) <- fs]
        ++ sieve ps q ((2*p,q):[(s,q-rem (q-y) s) | (s,y) <- fs])

A quick and dirty benchmarking function:

benchPrimes :: Num a => [a] -> IO ()
benchPrimes f = do start <- getCurrentTime
                   print . sum $ take 10000 f
                   print . flip diffUTCTime start =<< getCurrentTime

Time to run the benchmarks.

main :: IO ()
main = do benchPrimes primes_eratosthenes
          benchPrimes primes_euler

The sieve of Eratosthenes runs in about 0.02 seconds, while Euler’s sieve takes around 6.6 seconds, which means this implementation of Euler’s Sieve is about 300-400 times slower. Clearly you don’t want to use this version in practice.

Programming Praxis – Calculating Logarithms

December 18, 2009

In today’s Programming Praxis problem we have to implement Euler’s method for calculating logarithms (using Newton’s method for calculating square roots). Let’s get started.

First we define our desired accuracy.

eps :: Double
eps = 1e-7

Newton’s method is just repeating the same function until it’s close enough.

sqrt' :: Double -> Double
sqrt' n = head . filter (\x -> abs (x^2 - n) < eps) $
    iterate (\x -> x - (x^2 - n) / (2*x)) 1

In Euler’s algorithm we first find the smallest power of the base higher than n, and then keep calculating the geometric mean between the low and the high number, repeating for whichever of the two resulting intervals holds n until we’re close enough.

log' :: Double -> Double -> Double
log' b n = f 0 0 where
    f lo hi | b ** hi < n             = f lo (hi + 1)
            | abs (lo / hi - 1) < eps = arit
            | geom < n                = f arit hi
            | otherwise               = f lo arit
            where geom = sqrt' (b ** lo * b ** hi)
                  arit = (lo + hi) / 2

A quick test shows that everything is working correctly:

main = do print $ sqrt' 2
          print $ log' 10 612