Archive for April, 2010

Project Euler Problem #73!

April 26th, 2010

Problem #73 says:

Consider the fraction, n/d, where n and d are positive integers. If n

If we list the set of reduced proper fractions for d ≤ 8 in ascending order of size, we get:

1/8, 1/7, 1/6, 1/5, 1/4, 2/7, 1/3, 3/8, 2/5, 3/7, 1/2, 4/7, 3/5, 5/8, 2/3, 5/7, 3/4, 4/5, 5/6, 6/7, 7/8

It can be seen that there are 3 fractions between 1/3 and 1/2.

How many fractions lie between 1/3 and 1/2 in the sorted set of reduced proper fractions for d ≤ 12,000?

First off brute forcing this problem isn’t hard, and you can do it in well under a minute. However there is an elegant solution so readily available that it would be a shame not to go with it. The elegant method involves calculating a Farey Sequence (http://en.wikipedia.org/wiki/Farey_sequence) then counting the results. Alternatively we can just count how many results we would have without calculating the actual results.

Solution in f# and runs in ~.32 seconds. Requires .net 4 and fsharp.

#light
module ProjectEuler

open System
open System.Diagnostics

let watch = new Stopwatch()
watch.Start()

let rec FareySeq denom1 denom2 max =
    if (denom1+denom2) > max then
        0.0
    else
        1.0 + (FareySeq denom1 (denom1+denom2) max)+(FareySeq (denom1+denom2)denom2 max)

Console.WriteLine(FareySeq 3 2 12000)
watch.Stop()
Console.WriteLine(watch.Elapsed)
Console.Beep() |> ignore
Console.ReadKey() |> ignore

Tags: , ,
Posted in Project Euler | Comments (0)

Project Euler Problem #58!

April 25th, 2010

Problem #58 says:

Starting with 1 and spiraling anticlockwise in the following way, a square spiral with side length 7 is formed.

It is interesting to note that the odd squares lie along the bottom right diagonal,
but what is more interesting is that 8 out of the 13 numbers lying along both diagonals are prime;
that is, a ratio of 8/13 ≈ 62%.

If one complete new layer is wrapped around the spiral above, a square spiral with side length 9 will be formed.
If this process is continued, what is the side length of the square spiral for which the ratio of
primes along both diagonals first falls below 10%?

The below solution simply counts the offset between corners and checks that number to see if its prime or composite then the next cycle after its check it updates the offset amount and continues this way until its collected a list of prime and composite numbers that have less than 10% prime numbers. very simple, no optimizations applied. The solution could be made faster by using a simpler test and also by not keeping a list of composite numbers but a counter instead.

Solution in f# and requires the .net 4.0 run-time. runs in ~50 seconds.

#light
module ProjectEuler

open System
open System.Diagnostics
open Microsoft.FSharp.Linq
open Microsoft.FSharp.Collections
open Microsoft.FSharp.Math
open System.Numerics

let toBinary (n:BigInteger) =
    let mutable m = n
    let mutable r = []
    while m > BigInteger.Zero do
        r <- r @ [(m % (BigInteger 2))]
        m <- m / BigInteger 2
    r

let test (a:BigInteger) (n:BigInteger) =
    let (b:List<BigInteger>) = toBinary (n - BigInteger.One)
    let mutable d = BigInteger.One
    let mutable Prime = false
    let CheckList = [0 .. b.Length-1 ] |> List.rev
    for i in CheckList do
        let x = d
        d <- (d*d) % n
        if (d = BigInteger.One && x <> BigInteger.One && x <> n-BigInteger.One) then
            Prime <- true // complex
        if b.[i] = BigInteger.One then
            d <- (d*a) % n
    if d <> BigInteger.One then
        Prime <- true //complex
    Prime //if its still false then prime

let MillerRabin (n:uint64) (s:int) =
    let Rand = new System.Random()
    let mutable Prime = true
    if n < Convert.ToUInt64(Int32.MaxValue) then
        for j in [1 .. s+1] do
            let a = Rand.Next(1, Convert.ToInt32(n)-1)
            if (test (BigInteger a) (BigInteger n)) then
                Prime <- false
    else
        for j in [1 .. s+1] do
            let a = Rand.Next(1, Int32.MaxValue)
            if (test (BigInteger a) (BigInteger n)) then
                Prime <- false
    Prime

let watch = new Stopwatch()
watch.Start()

let rec Answer (currentnum:uint64) (step:uint64) (prime:List<uint64>) (composite:List<uint64>) (tick:bool) =
    let check = (float prime.Length)/((float composite.Length)+(float prime.Length))
    if check < 0.1 && currentnum > 100UL  then
        step-1UL
    else
        if (tick) then
            if (MillerRabin currentnum 10) then
                Answer (currentnum+step) step (currentnum::prime) composite (false)
            else
                Answer (currentnum+step) step (prime) (currentnum::composite) (false)
        else
            if (MillerRabin currentnum 10) then
                Answer (currentnum+step) (step+1UL) (currentnum::prime) composite (true)
            else
                Answer (currentnum+step) (step+1UL) (prime) (currentnum::composite) (true)

Console.WriteLine(Answer 3UL 2UL [] [1UL] true )
watch.Stop()
Console.WriteLine(watch.Elapsed)
Console.Beep() |> ignore
Console.ReadKey() |> ignore


Tags: , , ,
Posted in Project Euler | Comments (0)

Project Euler Problem #146!

April 25th, 2010

Problem #146 says:

The smallest positive integer n for which the numbers n^(2)+1, n^(2)+3, n^(2)+7, n^(2)+9, n^(2)+13, and n^(2)+27 are consecutive primes is 10.
The sum of all such integers n below one-million is 1242490.

What is the sum of all such integers n below 150 million?

This problem is rather complicated to say the least and my below solution is my first attempt at solving it, it takes about a half hour to run on a quad core machine. There are some tests in place to remove numbers as fast as possible for instance the numbers must be a multiple of 10 before we add a 1,3,7,9,13, or 27 to it, there is also a few modulus test being performed to narrow the search area down some. Then the list of possible numbers is filtered in parallel to make sure we only get consecutive primes. This solution uses the Miller-Rabin code posted earlier with a modified non-recursive toBinary method. There are 12 valid numbers.

Solution provided in f# and requires .net 4.0 and the Fsharp powerpack.

#light
module ProjectEuler

open System
open System.Diagnostics
open Microsoft.FSharp.Linq
open Microsoft.FSharp.Collections
open Microsoft.FSharp.Math
open System.Numerics

let toBinary (n:BigInteger) =
    let mutable m = n
    let mutable r = []
    while m > BigInteger.Zero do
        r <- r @ [(m % (BigInteger 2))]
        m <- m / BigInteger 2
    r

let test (a:BigInteger) (n:BigInteger) =
    let (b:List<BigInteger>) = toBinary (n - BigInteger.One)
    let mutable d = BigInteger.One
    let mutable Prime = false
    let CheckList = [0 .. b.Length-1 ] |> List.rev
    for i in CheckList do
        let x = d
        d <- (d*d) % n
        if (d = BigInteger.One && x <> BigInteger.One && x <> n-BigInteger.One) then
            Prime <- true // complex
        if b.[i] = BigInteger.One then
            d <- (d*a) % n
    if d <> BigInteger.One then
        Prime <- true //complex
    Prime //if its still false then prime

let MillerRabin (n:uint64) (s:int) =
    let Rand = new System.Random()
    let mutable Prime = true
    if n < Convert.ToUInt64(Int32.MaxValue) then
        for j in [1 .. s+1] do
            let a = Rand.Next(1, Convert.ToInt32(n)-1)
            if (test (BigInteger a) (BigInteger n)) then
                Prime <- false
    else
        for j in [1 .. s+1] do
            let a = Rand.Next(1, Int32.MaxValue)
            if (test (BigInteger a) (BigInteger n)) then
                Prime <- false
    Prime

let ntest n =
    let k = n / 10UL
    if (n%3UL=1UL || n%3UL=2UL)
        && n%13UL<>0UL
        && n%2UL=0UL
        && n%5UL=0UL &&
        ((n%7UL=3UL || n%7UL=4UL) &&
            (n%210UL=10UL ||
             n%210UL=80UL ||
             n%210UL=130UL ||
             n%210UL=200UL )) &&
        (k % 7UL = 1UL ||
         k % 7UL = 6UL ||
         k % 3UL <> 0UL ||
         k % 13UL <> 0UL ||
         k % 17UL <> 0UL ||
         k % 29UL <> 0UL ||
         k % 19UL <> 0UL ||
         k % 23UL <> 0UL ) then
        true
    else
        false

let nsquaretest n =
    let ns = n*n
    if ((MillerRabin ((ns)+13UL) 5) &&
        (MillerRabin ((ns)+3UL) 5) &&
        (MillerRabin ((ns)+7UL) 5) &&
        (MillerRabin ((ns)+9UL) 5) &&
        (MillerRabin ((ns)+1UL) 5) &&
        (MillerRabin ((ns)+27UL) 5) ) then
        if (MillerRabin (ns + 19UL) 5 ||
            MillerRabin (ns + 21UL) 5) then
            false
        else
            true
    else
        false

let rec sum (a:uint64) (n:List<uint64>)  =
    if n = List.empty then
        a
    else
        sum (a+n.Head) n.Tail

let watch = new Stopwatch()
watch.Start()

let check =
    [10UL .. 10UL .. 1000000UL]
    |> PSeq.filter ntest
    |> PSeq.toList

let Answer =
    check
    |> PSeq.filter nsquaretest
    |> PSeq.toList

Console.WriteLine(Answer |> sum 0UL)
Console.WriteLine(Answer.Length)
watch.Stop()
Console.WriteLine(watch.Elapsed)
Console.Beep() |> ignore
Console.ReadKey() |> ignore

Tags: , , ,
Posted in Project Euler | Comments (0)

FSharp Miller-Rabin primality test

April 20th, 2010

In preparation for some upcoming project euler problems I’ve been looking for faster methods of determining if a number is prime. One such method I’ve come across was the Miller-Rabin primality test which is a probabilistic test for determining if a given number is probably prime. The higher the value of s passed to the original function the higher the probability the answer is correct.

I found an example of this algorithm in python which I’ve rewritten to use in F#, I do not have the original python source but if you are interested in learning more you can check out the Wikipedia page.

The following code snippet uses .net 4′s BigInteger class for part of test. I have plans to rewrite it making more f# like right now its a direct translation of the python code.

open System.Numerics

let rec toBinary (n:BigInteger) r =
    if n > BigInteger.Zero then
        toBinary (n/(BigInteger 2)) (r @ [n% (BigInteger 2)])
    else
        r

let test (a:BigInteger) (n:BigInteger) =
    let (b:List<BigInteger>) = toBinary (n - BigInteger.One) []
    let mutable d = BigInteger.One
    let mutable Prime = false
    let CheckList = [0 .. b.Length-1 ] |> List.rev
    for i in CheckList do
        let x = d
        d <- (d*d) % n
        if (d = BigInteger.One && x <> BigInteger.One && x <> n-BigInteger.One) then
            Prime <- true // complex
        if b.[i] = BigInteger.One then
            d <- (d*a) % n
    if d <> BigInteger.One then
        Prime <- true //complex
    Prime //if its still false then prime

let Rand = new System.Random()

let MillerRabin (n:uint64) (s:int) =
    let mutable Prime = true
    for j in [1 .. s+1] do
        let a = Rand.Next(1, Convert.ToInt32(n)-1)
        if (test (BigInteger a) (BigInteger n)) then
            Prime <- false
    Prime

Tags: , ,
Posted in Math | Comments (0)