Sunday, November 21, 2010

Problem 50: Which prime, below one-million, can be written as the sum of the most consecutive primes?

Alright, we completed through problem 50!

First, I thought it was about time I implemented a decent prime sequence generator instead of continuing to combine Seq.initInfinite and isPrime. It was important to me that we create an infinite sequence of primes, rather than the typical finite prime sieve. I spent a considerable amount of effort playing with various algorithms and implementations of those algorithms. The final algorithm I settled on isn’t quite a sieve, but works on the general principle, allowing us to generate the millionth prime in about 4 seconds. Compare that to this infinite, true prime sieve which can generate the millionth prime in just over a second, but I prefer to stick with my algorithm, since I worked the core of it out independently. It may be worth pointing out that the memory used in my implementation is O(n), where n is the count of the prime, while any true sieve is O(p), where p is the prime itself.

To help us out, we finally get around to implementing an extension to the Seq module, infiniteRange (note that by marking the function inline, the compiler will automatically give a structural type signature to this function so that it can be used with int32, int64, bigint, float, etc.):

let inline infiniteRange start skip = 
    seq {
        let n = ref start
        while true do
            yield n.contents
            n.contents <- n.contents + skip

Next now we show our our primes implementation which is structurally generic as usual with integral type specific versions:

///p is "prime", s=p*p, c is "multiplier", m=c*p
type SievePrime<'a> = {mutable c:'a ; p:'a ; mutable m:'a ; s:'a}

///A cached, infinite sequence of primes
let inline primes_of (g:G<'a>) =
    let primeList = ResizeArray<_>()
    primeList.Add({c=g.three ; p=g.three ; m=g.three*g.three ; s=g.three*g.three})

    //test whether n is composite, if not add it to the primeList and return false
    let isComposite n = 
        let rec loop i = 
            let sp = primeList.[i]
            while sp.m < n do
                sp.c <-
                sp.m <- sp.c*sp.p

            if sp.m = n then true
            elif i = (primeList.Count-1) || sp.s > n then
                 primeList.Add({c=n ; p=n ; m=n*n ; s=n*n})
            else loop (i+1)
        loop 0

    seq { 
        yield g.two ; yield g.three

        for i in 1..primeList.Count-1 do
            yield primeList.[i].p

        yield! Seq.infiniteRange (primeList.[primeList.Count-1].p + g.two) g.two 
               |> Seq.filter (isComposite>>not)
let primes = primes_of gn
let primesL = primes_of gL
let primesI = primes_of gI

Now for problem 50. This was a pretty straight forward problem. We use an array in our algorithm but do not perform any mutation:

let problem50b =
    let primesArr = primes |> Seq.takeWhile ((>=)1000000) |> Seq.toArray
    let count p =
        let rec count p cur curMax =
            let start = primesArr.[cur]
            if start >= p/2 then curMax
                let rec consecutiveSum i sum =
                    if sum >= p then (i,sum)
                    else consecutiveSum (i+1) (sum + primesArr.[i])

                let i,sum = consecutiveSum cur 0            
                let curMax = max curMax (if sum = p then i-cur else 0)
                count p (cur+1) curMax
        count p 0 0

    |> Seq.maxBy count

Friday, November 19, 2010

Problem 49: What 12-digit number do you form by concatenating the three terms in this sequence?

This was a tough problem. Even after dreaming up an algorithm which would perform well enough, it went through several tweaks until I was satisfied. We resort to an array (though no mutation) and sequence expression in order to achieve short-circuiting and overcome F#’s lack of break and continue (again). This algorithm runs in about 1 minute and 500 milliseconds, just over the Project Euler suggested limit. I was able to uglify it a bit to get it to run in 59 seconds and 700 milliseconds, but I won’t even bother showing that here.

let problem49b =
    let fourDigitPrimes =
        |> Seq.filter isPrime
        |> Seq.toArray

    let arithmeticTriples = seq {
        for i in {0..fourDigitPrimes.Length-1} do
            for j in {i+1..fourDigitPrimes.Length-1} do
                for k in {j+1..fourDigitPrimes.Length-1} do
                    let a = fourDigitPrimes.[i]
                    let b = fourDigitPrimes.[j]
                    let c = fourDigitPrimes.[k]
                    //is arithmetic seq other than given in question
                    if b-a = c-b && a <> 1487 then 
                        yield (a,b,c)
    let arePerms (a,b,c) = 
        let setA = Set(Digits.fromInt a)
        let setB = Set(Digits.fromInt b)
        setA = setB && setB = Set(Digits.fromInt c)
    |> Seq.find arePerms
    |> (fun (a,b,c) ->
            [Digits.fromInt a; Digits.fromInt b; Digits.fromInt c] 
            |> Seq.concat 
            |> Digits.toInt64)

Saturday, November 13, 2010

Problem 48: Find the last ten digits of the series, 1^1 + 2^2 + 3^3 + ... + 1000^1000.

Again, easy problem using functional composition, the “hard” work was done in our previous implementation of the Digits module. The answer was produced surprisingly fast. I would like a better solution to finding the last 10 digits than to convert the entire sequence to an array.

let problem48a =
    let digits =
        |> Seq.sumBy (fun n -> BigInteger.Pow(BigInteger(n), n))
        |> Digits.fromBigInt
        |> Seq.toArray

    Array.sub digits (digits.Length - 10) 10
    |> Digits.toInt64 //just makes it easier to copy and paste result

Problem 47: Find the first four consecutive integers to have four distinct primes factors. What is the first of these numbers?

Not the best performing solution (perhaps inherently so, though), but a shining example of the power of functional composition (the hard work was done in implementing cfactorize in a previous problem, everything else is just composing built-in functional operators).

let problem47a =
    (Seq.initInfinite (fun i -> let n = i+3 in (n, cfactorize n)) //starting with 0+3=2
    |> Seq.windowed 4
    |> Seq.find (fun arr -> arr |> Seq.forall (fun (_,lst) -> lst |> List.length = 4))).[0]
    |> fst

Friday, November 12, 2010

Problem 46: What is the smallest odd composite that cannot be written as the sum of a prime and twice a square?

The difficult part of this problem was working around F#'s lack of break, continue, and return in for loops. hasForm would be most easily expressed using nested for loops but would require break, continue, and return to perform acceptably. Another desire is support for infinite range expressions, so that oddComposites could be written as {9..2..} |> Seq.filter (isPrime>>not).

let problem46b = 
    let form p s = p + 2*(s*s)

    //making sure to exclude 1    
    let oddComposites = 
        Seq.initInfinite id 
        |> Seq.filter (fun i -> i <> 1 && i%2=1 && isPrime i |> not)
    //cached for performance
    let primes = 
        Seq.initInfinite id 
        |> Seq.filter isPrime 
        |> Seq.cache

    let primesUpto n = 
        |> Seq.takeWhile ((>=) n)

    //think nested loops
    let hasForm n = 
        primesUpto n
        |> Seq.exists 
            (fun p ->
                |> (fun s -> form p s)
                |> Seq.find (fun r -> (r = n) || (r > n))) = n)

    |> Seq.find (hasForm>>not)