The essential problem is to determine the index i of the first Fibonacci number Fn to have at least d digits. The specific challenge is to determine the index of the first Fibonacci number to have 1000 digits.

Obviously, this problem could be brute forced by stepping through all the Fibonacci numbers and determining the answer. However, this is not very mathematically pleasing, nor is it very efficient. Instead, we will develop an algorithm to determine the answer using some nifty Fibonacci mathematics.

At this point, I'd like to introduce Mr. J. P. M. Binet. Binet was a French Mathematician and Physicist that lived and worked around the turn of the 19th century. Binet discovered a formula for calculating the nth Fibonacci number without having to know anything about the preceding Fibonacci numbers. His formula had been previously discovered by Euler and a fellow named Moivre, but the Formula is credited to Binet commonly. Binet's formula expressed the nth Fibonacci number in terms of the golden ratio (phi). phi = ( 1 + sqrt(5) ) / 2 and Fn( n ) = ( phi^n - ( -phi^-n ) ) / sqrt( 5 ) This is a tremendous time saving formula that allows instant calculation of any Fibonacci number. Interestingly enough, this formula always returns an integer even though it involves products of irrational numbers ( phi and sqrt(5) ). In any case, it should also be noted that if n is a real number, Fn will also be real number equal to some value between two Fibonacci numbers Fn( n ) < Fn( n + epsilon ) < Fn( n + 1 ) iff epsilon < 1 Another interesting property of this formula is that it can be reproduced perfectly with an approximation! Fn( n ) = round( phi^n / sqrt( 5 ) when n >= 0 This holds true because the inverse of phi is less than one, and its powers become very small ( and insignificant ) rapidly. So, in Binet's formula, the difference between phi^n and phi^-n is very nearly equal to phi when n becomes large. In fact, this difference is close enough to round correctly even when n is 0 or 1.

Now, since we have a simple formula for calculating Fn given n, we can simply solve for n to reverse the formula and determine the index of any fibonacci number.

Fn ~= phi^n / sqrt(5)
phi^n ~= Fn * sqrt(5)
log( phi^n ) ~= Fn * sqrt(5)
n log( phi ) ~= log( Fn * sqrt(5) )
n ~= log( Fn * sqrt(5) ) / log( phi )
[b]n ~= ( log( Fn ) + log( sqrt(5) ) ) / log( phi )[/b]

Now, if we round the result of this equation, we will always get the exact index of any Fibonacci number. However, there is some added value. In fact, if we use this formula on a non-Fibonacci number, it will give us an approximate Fibonacci index! n( Fn(n) ) < n( Fn(n) + eps ) < n( Fn(n+1) ) when eps < Fn(n+1) - Fn(n) If we then take the ceiling ( round up ) of the result of applying this formula to any number m, we can find the index of the first Fibonacci number greater than m. Lets call this function NextFib, and express it clearly: NextFib( m ) = ceiling( ( log( m ) + log( sqrt(5) ) ) / log( phi ) ) So, our final trick is to use this information to find the index of the first Fibonacci number with d digits. This is really quite trivial. We know that the first integer to have d digits is always: 10^(d-1) So, to find the index n of the first Fibonacci number Fn with d digits, we simply apply our magic formula: n( d ) = NextFib( 10^( d - 1 ) ) Now that we've developed a function to mathematically ( and elegantly, I might add ) find the index of the first Fibonacci number with d digits, we need to test it in code to verify that it indeed works. Here is some such code:

#include <cmath>
#include <iostream>

using namespace std;

inline double phi()
{
    return ( 1.0 + sqrt( 5.0 ) ) / 2.0;
}

/** Calculates the approximate Fibonacci number at index n */
double fib( double n )
{
    return pow( phi(), n ) / sqrt( 5.0 );
}

/** Calculates the exact integer Fibonacci number at index n */
int fibonacci( int n )
{
    return int( round( fib( (double)n ) ) );
}

/** Determines the approximate Fibonacci index for a given number ( Fibonacci or not ) */
double fibIdx( double Fn )
{
    return ( log(Fn) + log( sqrt(5.0) ) ) / log( phi() );
}

/** Determines the index of the next Fibnacci number to have d digits */
int magic( int d )
{
    return int( ceil( fibIdx( pow( 10.0, d-1.0 ) ) ) );
}


int main(int argc, char *argv[])
{
    cout << "Determining indices (i) of first Fibonacci numbers (Fn) to have some number of digits (d)" << endl;
    cout << "=========================================================================================" << endl;
    cout << "d in [2:10)" << endl;
    cout << "-----------" << endl;
    int i, d;
    for( d=2; d<10; d++ )
    {
        i = magic( d );
        cout << "  d=" << d << " i=" <<  i << " Fn(i)=" << fibonacci(i) << endl;
        cout << "      Fn(i-1)=" << fibonacci( i-1 ) << "  Fn(i+1)=" << fibonacci( i+1 ) << endl;
    }

    cout << "----------------------------" << endl;
    cout << "d in [10:100) by steps of 10" << endl;
    cout << "----------------------------" << endl;
    for( d=10; d<100; d+=10 )
    {
        i = magic( d );
        cout << "  d=" << d << " i=" <<  i << " Fn(i)=" << fib(i) << endl;
        cout << "      Fn(i-1)=" << fib( i-1 ) << "  Fn(i+1)=" << fib( i+1 ) << endl;
    }

    cout << "------" << endl;
    cout << "d=1000" << endl;
    cout << "------" << endl;
    d = 1000;
    i = magic( 1000.0 );
    if( i < 0 )
        cout << "overflow occurred" << endl;
    else
    {
        cout << "  d=" << d << " i=" <<  i << " Fn(i)=" << fib(i) << endl;
        cout << "      Fn(i-1)=" << fib( i-1 ) << "  Fn(i+1)=" << fib( i+1 ) << endl;
    }
    cout << "=========================================================================================" << endl;

    return 0;
}

And the results:

Determining indices (i) of first Fibonacci numbers (Fn) to have some number of digits (d)
=========================================================================================
d in [2:10)
-----------
  d=2 i=7 Fn(i)=13
      Fn(i-1)=8  Fn(i+1)=21
  d=3 i=12 Fn(i)=144
      Fn(i-1)=89  Fn(i+1)=233
  d=4 i=17 Fn(i)=1597
      Fn(i-1)=987  Fn(i+1)=2584
  d=5 i=21 Fn(i)=10946
      Fn(i-1)=6765  Fn(i+1)=17711
  d=6 i=26 Fn(i)=121393
      Fn(i-1)=75025  Fn(i+1)=196418
  d=7 i=31 Fn(i)=1346269
      Fn(i-1)=832040  Fn(i+1)=2178309
  d=8 i=36 Fn(i)=14930352
      Fn(i-1)=9227465  Fn(i+1)=24157817
  d=9 i=40 Fn(i)=102334155
      Fn(i-1)=63245986  Fn(i+1)=165580141
----------------------------
d in [10:100) by steps of 10
----------------------------
  d=10 i=45 Fn(i)=1.1349e+09
      Fn(i-1)=7.01409e+08  Fn(i+1)=1.83631e+09
  d=20 i=93 Fn(i)=1.22002e+19
      Fn(i-1)=7.54011e+18  Fn(i+1)=1.97403e+19
  d=30 i=141 Fn(i)=1.31151e+29
      Fn(i-1)=8.10559e+28  Fn(i+1)=2.12207e+29
  d=40 i=189 Fn(i)=1.40987e+39
      Fn(i-1)=8.71347e+38  Fn(i+1)=2.28122e+39
  d=50 i=237 Fn(i)=1.5156e+49
      Fn(i-1)=9.36695e+48  Fn(i+1)=2.4523e+49
  d=60 i=284 Fn(i)=1.00694e+59
      Fn(i-1)=6.22325e+58  Fn(i+1)=1.62927e+59
  d=70 i=332 Fn(i)=1.08246e+69
      Fn(i-1)=6.68997e+68  Fn(i+1)=1.75146e+69
  d=80 i=380 Fn(i)=1.16364e+79
      Fn(i-1)=7.19168e+78  Fn(i+1)=1.88281e+79
  d=90 i=428 Fn(i)=1.25091e+89
      Fn(i-1)=7.73103e+88  Fn(i+1)=2.02401e+89
------
d=1000
------
overflow occurred
=========================================================================================

We can see that for d = [2:10) , [10:100) the calculations were correct. It indeed found many indices for Fibonacci numbers with d digits. Sadly, the math overflowed when trying to calculate the numbers with d=1000.

Rather than trying to find a way to make this work in c++, I skipped over to python ( which sometimes has better numerical resolution ) and tried it. Fortunately, Python was able to calculate log_10( 10*999 ) while c++ was not.

from math import *

phi = ( 1 + 5**.5 ) / 2

def fib( n ):
	return phi**n / 5**.5

def idx( Fn ):
	return ( log( Fn ) + log( 5**.5 ) ) / log( phi )

def magic( d ):
	return int( ceil( idx( 10**(d-1) ) ) )

i=magic( 1000 )
print i

And, python reported that the index was 4782, though the fibonacci number for this index couldn't be calculated using Binet's formula due to an overflow:

d=1000, i(d)=4782
Traceback (most recent call last):
  File "tester.py", line 16, in <module>
    print "Fn(%i)=%d" %( i, fib(i) )
  File "tester.py", line 6, in fib
    return phi**n / 5**.5
OverflowError: (34, 'Numerical result out of range')

So, there it is. The maths for this solution are not particularly complicated, nor is the logic, nor is the code. However, I think that the results are quite pleasing aesthetically and mathematically.

Close but you still did not need to program for this problem.

We have this wonderful formula :

Fib(n) = round( Phi^n/√5)


We can manipulate the above equation to get the number of digit we are looking for by using logarithm.


Fib(n) = log( phi^n / phi(√5) )

Lets denote Fib(n) to be X, the number of digit we are looking for,
then we have :

x = log( phi^n / phi(√5) )

x = log(phi^n) - log(√5)

x = n*log(phi) - 0.5*log(5)

( x + 0.5*log(5) ) / log(phi) = n

Thus we have a formula for n :

n = ( x + 0.5*log(5) ) / log(phi)

Now we want to digit to be 1000. But in the original equation
we have the round function. So to cancel that out, we subtract
1 from the input x and ceil our final equation.
Thus our final equation is :

n = ceil ( ( (x - 1) + 0.5*log(5) ) / log(phi) )

Now plugging x = 1000 and phi = 1.618, we get :

n = 4782.

Ah, yes.

In developing the system of equations, I overlooked one simple reduction

Given:

Fn( m ) = ceiling( ( log( m ) + log( sqrt(5) ) ) / log( phi ) )
n( d ) = Fn( 10^( d - 1 ) )

If we select a base 10 logarithm for the Fn function, and collapse the two functions together we get: n( d ) = ceiling( ( log_10( 10^(d-1) ) + log_10( sqrt(5) ) ) / log_10( phi ) ) Since we know that log_10( 10^k) == k, we can simply reduce the function to: n( d ) = ceiling( ( d - 1 + log_10( sqrt(5) ) ) / log_10( phi ) ) Finally, if we really don't like finding the square root of 5 every time, we can simplify that logarithm by log_b( a^x ) == x * log_b( a ) to get: n( d ) = ceiling( ( d - 1 + 0.5 * log_10( 5 ) ) / log_10( phi ) ) Of course, this equation should be robust against overflow for VERY large numbers indeed. It is also the equation you arrived at.

This article has been dead for over six months. Start a new discussion instead.