Problem #66 says:

Consider quadratic Diophantine equations of the form:

x^(2) – Dy^(2) = 1

For example, when D=13, the minimal solution in x is 649^(2) – 13×180^(2) = 1.

It can be assumed that there are no solutions in positive integers when D is square.

By finding minimal solutions in x for D = {2, 3, 5, 6, 7}, we obtain the following:

3^(2) – 2×2^(2) = 1
2^(2) – 3×1^(2) = 1
9^(2) – 5×4^(2) = 1
5^(2) – 6×2^(2) = 1
8^(2) – 7×3^(2) = 1

Hence, by considering minimal solutions in x for D ≤ 7, the largest x is obtained when D=5.

Find the value of D ≤ 1000 in minimal solutions of x for which the largest value of x is obtained.

So there are several ways to tackle this problem at first glance; Just role a pair of nested loops and try all the possible combinations until you get a working one then keeping track of the max value for X and D, the second is to find something else because the problems after number 10 are not conducive to a brute force approach. The first option isn’t feasible because the number of X,Y pairs you have to check to find a solution is exponentially greater than it appears at first glance. So we need to find something else.

What we can find with a simple Wikipedia search is that this particular form of equation is called a “Pell Equation” (Oddly enough it seems to be misnamed — look it up). This gives us some information that was already in the problem statement and a method for attacking the problem more efficiently. We’ve read in two different sources now that we don’t want D to be a perfect square; so the first thing we need is a simple test to check if a number is a perfect square.

//This is separated out because it's useful in other problems too.
let PerfectSquare (d:int) =
    let sqrtd = Math.Sqrt((double d))
    if ((double d) %  sqrtd = 0.0) then
        true
    else
        false

Just grab the square root of a number and then take the number mod it’s square root. if you end up with no remainder you have a perfect square and we want to skip that number. There might be faster methods but for this solution and the amount of time we call it this works well enough. So what next? Well the solution I went with is the continued fraction method; so we’ll need a method to generate continued fractions. Specifically we want to generate a Continued fraction set going towards the square root of the D value we’re currently checking. Then those values are used to calculate the fundamental solution to the Pell Equation which gives us the minimal X value for a given D value.

Full solution provided in F# and was developed with the F# mono-develop plug-in. Runs in under a second and prints the D value and the maximum X value; the X value should explain why using a nested loop is a really bad idea.

open System
open System.Numerics
open System.Diagnostics
open System.Collections.Generic

let GenerateCFraction (target:int64) =
    let e = new System.Collections.Generic.List<int64>()
    let SqrtTarget = int64 (Math.Sqrt(double target))
    let a = ref SqrtTarget
    let p = ref 0L
    let q = ref 1L
    e.Add(!a)
    p := !a * !q - !p
    q := (target - !p * !p) / !q
    a := (!p + SqrtTarget) / !q
    while !q <> 1L do
        e.Add(!a)
        p := !a * !q - !p
        q := (target - !p * !p) / !q
        a := (!p + SqrtTarget) / !q
    e.Add(!a)
    e

let SolvePell (cFrac:System.Collections.Generic.List<int64>) =
    let l = cFrac.Count - 1
    let per =
        if l % 2 = 0 then
            l - 1
        else
            2 * l - 1
    let a  = ref (BigInteger cFrac.[0])
    let a1 = ref 1I
    let b  = ref 1I
    let b1 = ref 0I
    let t = ref 0I
    a := (BigInteger cFrac.[0])
    b  := ((BigInteger cFrac.[0]) * !b + !b1)
    for i = 1 to (per) do
        t := !a
        a := (BigInteger cFrac.[(i-1)% l + 1]) * !a + !a1
        a1 := !t
        t := !b
        b := (BigInteger cFrac.[(i-1)% l + 1]) * !b + !b1
    !a

let PerfectSquare (d:int) =
    let sqrtd = Math.Sqrt((double d))
    if ((double d) %  sqrtd = 0.0) then
        true
    else
        false

[<EntryPoint>]
let main args =
    let sw = new Stopwatch()
    sw.Start()
    let x = ref 0I
    let maxx = ref 0I
    let maxd = ref 0
    for d = 2 to 1001 do
        if (PerfectSquare d) then
            ignore
        else
            let c = GenerateCFraction (int64 d)
            x := (SolvePell c)
            if (!x > !maxx) then
                maxx := !x
                maxd := d
            ignore
    sw.Stop()
    Console.WriteLine sw.Elapsed
    Console.WriteLine !maxd
    Console.WriteLine !maxx
    0

Problem #87 says:

The smallest number expressible as the sum of a prime square, prime cube, and prime fourth power is 28. In fact, there are exactly four numbers below fifty that can be expressed in such a way:

28 = 2^(2) + 2^(3) + 2^(4)
33 = 3^(2) + 2^(3) + 2^(4)
49 = 5^(2) + 2^(3) + 2^(4)
47 = 2^(2) + 3^(3) + 2^(4)

How many numbers below fifty million can be expressed as the sum of a prime square, prime cube, and prime fourth power?

Not a whole lot of complicated thought here. Just generate a list of primes from 2 to the square root of (50000000 – 24) and nest some loops. Simple enough, the only other thing we need to do is remove the duplicates in the resulting list and count the entries. It’s a rather easy problem considering where it falls on the list of problems sorted by difficulty on the project euler website. It might be faster to zip or collect the lists in an intermediate form instead of just nested loops; but this solution is fast enough for me with a run time of 16 seconds.

In other news I’m really starting to get into this loop & yield pattern. It’s a really concise way of generating a set of results, makes it easy to control the element types. Plus it seems a lot more ‘functional’ than my previous system. Maybe I’m just getting used to seeing it when I read f# code so I’ve started using it myself; sorta like learning through osmosis…

Anyways, solution provided in F# written in monodevelop. Also this makes the 90th problem I’ve solved on project euler only 225 more!

open System
open System.Collections.Generic
open System.Diagnostics

let rec sieve = function
    | (p::xs) -> p :: sieve [ for x in xs do if x % p > 0L then yield x ]
    | []      -> []

let KnownPrimes = sieve [2L .. 7071L]

let Numbers max =
    List[
        for x in KnownPrimes do
            for y in KnownPrimes do
                for z in KnownPrimes do
                    let candidate = (x*x) + (y*y*y) + (z*z*z*z)
                    if candidate < max then
                        yield candidate
    ]
    |> Seq.distinct
    |> Seq.toList

[<EntryPoint>]
let main args =
    let sw = new Stopwatch()
    sw.Start()
    let ans = Numbers (int64 args.[0])
    printfn "%i" ans.Length
    sw.Stop()
    Console.WriteLine(sw.Elapsed)
    0

Problem #85 says:

By counting carefully it can be seen that a rectangular grid measuring 3 by 2 contains eighteen rectangles; Although there exists no rectangular grid that contains exactly two million rectangles, find the area of the grid with the nearest solution.

The project euler site has a nifty picture showing an example. This solution is really simple; in fact I would be willing to bet that anybody who has taken high school algebra should be able to solve this with an efficient formula. The formula we’ll use is based off of the formula for adding the numbers from one to x which looks like this ((x*x)+x)/2 so what do we do when we have two dimensions? Simply multiply to obtain (((x*x)+x)/2)(((y*y)+y)/2) which becomes (((x*x)+x)((y*y)+y))/4 since this grows pretty fast we can cap the endpoints of our search to 100 for both x and y. Another thing to note is the below solution will give us the same rectangle rotated 90 degrees which I simply ignore.

Solution provided in f# and runs in under a second. Written using mono and the f# mono-develop plug-in!

open System
open System.Collections
open System.Diagnostics
open Microsoft.FSharp.Collections

let RectInRect target =
    seq{
    for x = 1 to 100 do
        for y = 1 to 100 do
            let rectangles  = (((x*x)+x)/2)*(((y*y)+y)/2)
            if (float rectangles) > (float target) * 0.9995 && (float rectangles) < (float target) * 1.0005 then
                yield (x,y)
    }
    |> Seq.toList

[<EntryPoint>]
let main args =
    let sw = new Stopwatch()
    sw.Start()
    let r = RectInRect 2000000
    for x in r do
        Console.WriteLine ((fst x) * (snd x))
        Console.WriteLine (Convert.ToString(fst x)+ " " + Convert.ToString(snd x))
    sw.Stop()
    Console.WriteLine(sw.Elapsed)
    0

Problem #76 says:

It is possible to write five as a sum in exactly six different ways:

4 + 1
3 + 2
3 + 1 + 1
2 + 2 + 1
2 + 1 + 1 + 1
1 + 1 + 1 + 1 + 1

How many different ways can one hundred be written as a sum of at least two positive integers?

What we’re looking at here is called integer partitioning. A simple process where one counts (or lists) the various ways one can add numbers to hit a target number. My first solutions to this problem was to write a recursive algorithm that generated all possible combinations then just count the resulting list. It worked fine for the example of 5 and for the test values of 10, 15, and 20 but each time it doubled the run time. So anybody who’s ever heard of the chessboard problem knows that this solution is no good because it wouldn’t finish in my lifetime. So I did some searching and found an algorithm for generating the solution using an iterative process and was able to get it solved in no time at all.

I wrote this in f# using mono (and the new mono-develop plug-in for f#!!!) solves the problem in about 6 seconds. I’ve left the code that returns a list of the partitions commented out in case I have need of it later but be careful if you re-enable it as it consumes memory very fast (it hits 10GB of ram usage before solving the problem) as it stands right now it just increments a counter when it would output a partition.

Another thing that is making a new appearance with this solution is the usage of reference cells, something I haven’t used extensively before on this site; however they were required to get the algorithm because you can’t use mutable variables in all instances for more information check the msdn site.

#light
open System
open System.Linq
open System.Collections.Generic
open System.Diagnostics

let Partitions n =
    //List[
    let ans = ref 0
    let a = new List<int>(n+1)
    a.AddRange(List.init (n+1) (fun x -> 0))
    let k = ref 1
    a.[0] <- 0
    let y = ref (n-1)
    while !k<>0 do
        let x = ref (a.[!k-1] + 1)
        k := !k-1
        while (2 * !x <= !y) do
            a.[!k] <- !x
            y := !y - !x
            k := !k+1
        let l = ref(!k+1)
        while (!x <= !y) do
            a.[!k] <- !x
            a.[!l] <- !y
            ans:=!ans+1//yield a.Take(!k+2).ToList
            x := !x +  1
            y := !y-1
        a.[!k] <- !x + !y
        y := !y + !x - 1
        ans:=!ans+1//yield a.Take(!k+1).ToList
    !ans
    //]

[<EntryPoint>]
let main args =
    let sw = new Stopwatch()
    sw.Start()

    let n = (Partitions (System.Int32.Parse args.[0])) - 1
    Console.WriteLine (Convert.ToString n)
    sw.Stop()
    Console.WriteLine(sw.Elapsed)
    0

Now that you’ve seen that solution I’ll give you a second solution. You can solve this problem without any programing if you head over to wolfram alpha and get to know the PartitionsP command and can figure out how many unneeded solutions you have to remove from its result to get the correct answer (Hint: problem says sum of two numbers!)

I had some free time between classes the other day and decided to take a look at Haskell. Mostly because I’ve been hearing good things about the language and because I’ve been learning f# and looking at a second functional language can’t hurt right? And I’ve decided that the people that wrote Haskell did it Specifically to compete in code golf competitions. It’s a really compact language in my opinion, maybe too compact. But I looked at several solved project euler problems in Haskell and tried to create the simplest solution I could think of. All of these were created in ghci.

So my favorite feature of Haskell (so far I’m still learning obviously) list comprehension makes this problem laughably simple basically we tell the language what we want and really nothing else. there is nothing in this code snippet (which is a complete solution) that doesn’t contribute to the answer. This also holds for problem #6, again this is remarkably simple using list comprehension. The map function gives us a simple way to perform an action on a set, and when you combine map with list comprehension you can do some reasonably complex stuff in a really simple manner.

Anyway onto the solutions. All solutions can be run in ghci and all run under the one minute mark. Most of these problems show up in other languages elsewhere on this site if your curious.

//Project euler problem #1
sum [x | x <- [1..1000-1], x `mod` 3 == 0 || x `mod` 5 == 0]

//Project euler problem #6
((sum [ x | x <- [1..100]]) ^2)-(sum [ x*x | x <- [1..100]])

//Project euler problem #16
import Data.Char
let twotothousand = show (2^1000)
let twoetcchars = map digitToInt twotothousand
sum twoetcchars

//Project euler problem #20
let factorial n = product [ 1 .. n]
let digits = show (factorial 100)
let sumdigits = sum (map digitToInt digits)
sumdigits

//Project euler problem #48
(sum [x^x | x<- [1..1000]] ) `mod` 10000000000

//Project euler problem #56
let checks = [a^b | a<- [1..100] , b <- [1..100]]
let digitsum n = sum (map digitToInt (show n))
let finalcheck = [digitsum x | x <- checks]
maximum finalcheck

Though at this point I have to say, I prefer F# to Haskell. Probably because I’ve got more experience with F# but Haskell is interesting enough I might keep looking into it. I will say that I wish it was as easy to create an infinite list in F# as it is in Haskell. In Haskell one can do this accidentally but in F# its got a weird syntax. Anyways thats what I got for now.