Project Euler 57: Investigate the expansion of the continued fraction for the square root of two.

Project Euler 57: Investigate the expansion of the continued fraction for the square root of two.

Written by on 20 August 2011

Topics: Project Euler

I found Problem 57 of Project Euler to be a rather interesting problem, with more than one solution. The problem description reads

It is possible to show that the square root of two can be expressed as an infinite continued fraction.

√ 2 = 1 + 1/(2 + 1/(2 + 1/(2 + … ))) = 1.414213…

By expanding this for the first four iterations, we get:

1 + 1/2 = 3/2 = 1.5
1 + 1/(2 + 1/2) = 7/5 = 1.4
1 + 1/(2 + 1/(2 + 1/2)) = 17/12 = 1.41666…
1 + 1/(2 + 1/(2 + 1/(2 + 1/2))) = 41/29 = 1.41379…

The next three expansions are 99/70, 239/169, and 577/408, but the eighth expansion, 1393/985, is the first example where the number of digits in the numerator exceeds the number of digits in the denominator.

In the first one-thousand expansions, how many fractions contain a numerator with more digits than denominator?

First thing I will present is a brute force solution, which for all practical purposes are fast enough. The second solution is a closed form approximation to the problem, so it can be solved as fast as you can punch the calculator.

Brute Force Solution

The approach to the brute force solution is to do a bit of pattern recognition of the fractions, a bit like the once you get on an IQ test, please continue the this sequence….

We have

From this we can derive the following two formulas

 

 

where n means the numerator, and d means the denominator.  To calculate this in C# we would need two temporary variables to store the old values of the numerator and denominator, but it can be rewritten to

which means we only need two variables in order to store it. For all practical purposes this means nothing for the performance, but it makes the code a bit cleaner.
To check the number of digits we can use the log10 function and compare the integer parts. The total code for this looks like

const int limit = 1000;
int result = 0;

BigInteger den = 2;
BigInteger num = 3;

for (int i = 1; i < limit; i++) {
     num += 2 * den;
     den = num - den;
     if ((int)BigInteger.Log10(num) >(int)BigInteger.Log10(den))
         result++;
}

Unfortunately the expansions result in some rather huge numbers, and thus we have to use BigInteger if we want to get the result. This slows down the performance a bit which shows from the 8ms the solution took:

The number of iterations with more digits in numerator is 153
Solution took 8 ms

Approximate Solution

We know from the problem description that the first time we encounter a numerator with more digits than the denominator is after 8 expansions. For some reason this got me a bit curious to see when the next solution would arise. This happens after expansion 13 and then after 21 and 26. So after a few iterations a pattern started emerging that we would have 8,5,8,5,8,5 expansions between the solutions to the problem. I tried to iterate several more times and got the following

1234567x1234x1234567x1234x
1234567x1234x1234567x
1234567x1234x1234567x1234x
1234567x1234x12x1234x
1234567x1234x1234567x1234x
1234567x1234567x1234x
1234567x1234x1234567x1234x
1234567x1234567x1234x
1234567x1234x1234567x
1234567x1234x1234567x1234x
1234567x1234x1234567x

which shows that there is a trend of 13 solutions per 13 expansions, but the pattern is not completely regular. However, for fun and giggles lets try to assume that the pattern was regular and see how close we get. We could form a solution looking like

where the x is the limit and the means floor.

Plotting the error term of this approximation against the actual value of problem shows that the approximation is within +/- 1 of the actual solution.

In 68% of the cases between 1 and 1000 it hits the target, in 8,4% it shoots 1 above and in 23,6% it is just one short. For the value 1000 it hits spot on. So even if it is not a completely accurate approximation it is rather close.

Wrapping up

I presented you with two solutions, an accurate and rather fast brute force approach and an approximate but closed form solution. It might be possible to enhance the closed form solution, I haven’t investigated that in depth yet.

As usual you can find the source code, if you want a close peek at it.

How did you approach the problem?

12 Comments For This Post I'd Love to Hear Yours!

  1. SuprDewd says:

    I found the same formula for generating n/d using the last n/d after inspecting the problem for a bit.
    And then I just looped 1000 times and incremented when log10(n) > log10(d).

    Exactly the same as you did. 🙂

  2. Khaur says:

    I found different (but equivalent) formulas for N and D: . Now for D this is a second order linear recurrence, which means we can get closed forms:

    However I didn’t manage to find a clever way to use those closed forms to solve the problem at hand (even for brute forcing they are less efficient than the recurrence formulas).

  3. Kristian says:

    I don’t mean to censor you, but I deleted your second comment and updated the second formula so it parses.

    I have the exact same issues with the latex formulas as you. I just have the ability to edit 🙂

  4. Kristian says:

    Nicely found. This would be a good thing if you had to answer for one specific number. But you are right, since we need to check every number anyway, closed form or recurrence does not make much of a difference.

  5. sash_shea says:

    Hello there i found this amazing pattern.
    The first 8 sequence is S[3/2, 7/5, 17/12, 41/29, 99/70, 239/169, 577/408, and 1393/985].
    S2 is formed from S1, S1 numerator + S1 denominator is S2’s denominator.
    S2’s denominator plus S1’s denominator is S2’s numerator.

  6. S.Ged says:

    Very nice approximate !!!
    Thanks

  7. Akshat says:

    I have written the code in C++. The answer it is giving me is 123 which is wrong. Can anyone help me with it?

    #include
    #include
    using namespace std;

    int main()
    {
    long double d1=2,n1=3;
    int count=0;
    for (int i=0;i(long)log10(d1))
    {
    count++;
    cout<<"num is "<<n1<<" and den is "<<d1<<"\n";
    }
    }
    cout<<count;
    }

  8. Lou Weed says:

    Pretty sure your brute force solution will get a significant speed up if you do this instead of performing a log on 2000 numbers (some of which are ridiculously big):


    proc solve {{n 1000}} {
    puts "57. In the first $n expansions of the square root of 2, how many fractions contain a numerator with more digits than denominator?"

    set answer 0
    set num 3
    set den 2
    set power10 10

    for {set expansion 1} {$expansion=$power10} {
    # if the denominator is still less than 10/100/1000/etc. then increase the answer
    if {$den<$power10} { incr answer }
    set power10 [expr {$power10*10}]
    }

    # next expansion
    set num [expr {$num+2*$den}]
    set den [expr {$num-$den}]
    }

    return $answer
    }

    So most iterations, the only thing I do is check to see if the numerator has gone over the next power of 10, which is a simple comparison. On the occasions when it has, check if the denominator made it over too. If not, add to the answer.

  9. Lou Weed says:

    That code’s Tcl btw

  10. Jean-Marie Hachey says:

    Table 1
    Number of fractions containing a numerator with more digits
    than the denominator (NFND)
    Re: Project Euler Problem 57

    http://img11.hostingpics.net/pics/838263pe57tab1.jpg

    Table 1 shows that in the first 1000 expansions, 153 fractions contain
    a numerator with more digits than the denominator; this quantity reaches a value of 1508 fractions after 10 000 expansions.
    Slopes of consecutive points also show that the points are not collinear.

    ___

    Table 2
    Distribution of digits in the numerator and denominator.
    Re: Project Euler Problem 57

    http://img11.hostingpics.net/pics/138849pe57tabletwo.jpg

    In Table 2, as predicted with the algo, we find 6 fractions (bold) where the number of digits in the numerator is greater.

    ___

    Sources:

    1) Project Euler 57
    Kristian’s algorithm
    http://www.mathblog.dk/files/euler/Problem57.cs
    2) Microsoft Visual C# 2010 Express
    (Reference added : System.Numerics)
    3) Microsoft Office Excel (2007)
    4) General interest:
    The Square Root of Two to 1 Million Digits
    http://apod.nasa.gov/htmltest/gifcity/sqrt2.1mil

  11. Rahul says:

    Hey Kristian,
    The pattern you are looking for is Parabola. The repeating sequence 8,5,8,5… will make a parabola (Or hyperbola) like shape.

    Graph

  12. Rahul says:

    Last update was my bad
    This is the sequence (tested upto 10K)

    (8,5,8,5,8,5,8,8,5,8,5,8,5,3)
    (5,8,5,8,5,8,8,5,8,5,8,5,8,8,5,8,5,8,8,5,8,5,8,5,8,8,5,8,5,8,5,3)*x

    Graph

Leave a Comment Here's Your Chance to Be Heard!

You can use these html tags: <a href="" title=""> <abbr title=""> <acronym title=""> <b> <blockquote cite=""> <cite> <code> <del datetime=""> <em> <i> <q cite=""> <s> <strike> <strong>

You can use short tags like: [code language="language"][/code]

The code tag supports the following languages: as3, bash, coldfusion, c, cpp, csharp, css, delphi, diff,erlang, groovy, javascript, java, javafx, perl, php, plain, powershell, python, ruby, scala, sql, tex, vb, xml

You can use [latex]your formula[/latex] if you want to format something using latex

This site uses cookies.