4 Ways to Generate the Fibonacci Numbers

In keeping with our recent theme of generating famous number structures, we'll consider the Fibonacci sequence today. This post was also inspired by some of the examples and exercises in the Structure and Interpretation of Computer Programs (namely Exercise 1.13 and Exercise 1.19 for those that came for the Scheme solutions.)

The story goes that the Fibonacci numbers were so named by the Italian mathematician Leonardo Fibonacci in the 1200's, when he discovered the sequence while forming a mathematical model for rabbit populations. However evidence suggests that the sequence was known to others before then (it seems to appear in earlier Indian mathematics, for instance.)

The Fibonacci numbers also have an interesting relationship with the Golden Ratio, namely that as we take larger and larger numbers in the sequence, the ratio between one number and the next approaches the Golden Ratio.

Or mathematically:

And these numbers pop up in all sorts of places. Remember Pascal's triangle? Turns out if you draw some “shallow diagonal” lines across the triangle and sum up the numbers on each line, you get the Fibonacci sequence too!

(There's a nice picture of this here.)

Though that history is all well and good, computation is always much more interesting! We'll look at four different ways to generate the Fibonacci numbers, with Clojure taking centre stage as our JVM language today, including some guest appearances from Python (or Jython if you like, so as to keep with the JVM theme.)

So without further ado...

1. Using a Recursive Process

The simplest way to generate the Fibonacci numbers is recursively.

Let n be some positive integer. Then the recurrence relation, including the base cases, for the nth Fibonacci number is as follows:

We can express this exactly in Clojure code with a recursive function:

We can then call this function from the Clojure REPL to generate, say, the first 20 non-zero Fibonacci numbers:

     (map fibonacci (range 1 (inc 20)))
     ;> (1 1 2 3 5 8 13 21 34 55 89 144 233 377 610 987 1597
        2584 4181 6765)

This is a concise and elegant way to compute the Fibonacci numbers. However, using recursion in this way is actually so inefficient that it was worthy of a quip in Code Complete!

To understand why, notice that we've got both of the following:

• overlapping subproblems, ie. When we compute F₅ for example, we compute F₄ and F₃ , and F₄ itself computes F₃ again in its own branch. Both of these calls to F₃ make their own separate calls to F₂ and F₁ , and you can see how this quickly gets ridiculous. So we're doubling (and tripling, and quadrupling, …) up on the work we end up doing for each calculation!
• delayed operations. This means that the code cannot possibly be optimised to be tail-recursive, so we're creating a new stack frame for each of those calls!

The inefficiency of this method was also explored in SICP, of course. The process has order of growth O(2ⁿ) in time and O(n) in space.

(Try computing (fibonacci 500) if you like! Funny story: I was originally going to do this in the hopes that it would quickly give me a StackOverflow exception that I could talk about here for a bit. But an hour later the thing still hadn't finished and I gave up.)

So though the recursive method has a certain beauty to it, it's probably not a good way to calculate Fibonacci numbers. Luckily, there are more efficient means of computation!

2. Using an Iterative Process

A better approach is to compute the Fibonacci numbers using an iterative process. Intuitively, this means that we build up the next numbers at each step (instead of branching down), which means we don't have any repeated calculations this time.

Basically, at each point we're applying a transformation to a pair of integers [current, next]:

• we start off with the initial values current = 0, next = 1
• at each point, we apply the transformation 
  [current, next] → [next, current+next].

So for example, [0,1] → [1,1] → [1,2] → [2,3]. 

After n transformations, the nth Fibonacci number appears in the current (left-hand) slot.

Now most programming languages these days have iteration front and centre, and so are well suited to representing these particular 'looping variable update' constructs. In the same way that the recursion method translated easily into Clojure, the iterative method translates very easily into any of these languages.

Here's the iterative method in Python (or Jython, its JVM variant):

Generating the first 20 non-zero Fibonacci numbers in the Python interpreter:

     map(fibonacci, range(1,21))
     #> [1, 1, 2, 3, 5, 8, 13, 21, 34, 55, 89, 144, 233, 377,

        610, 987, 1597, 2584, 4181, 6765]

The numbers agree with our recursive definition from earlier, so that's good.

In Clojure, we can use the loop-recur construct to generate an iterative process by using the aforementioned tail-recursion. Syntactically it looks very similar to how it would look in a procedural language, but behind the scenes we are actually using recursion to run the iteration, the state transformations happen as the recursive calls change the values of the input parameters each time. Since the recursive call is in the 'tail' position (that is, the call is the last thing that happens and isn't nested inside another expression), the Java Virtual Machine can optimise the recursion away, so we don't accumulate stack frames as we go.

Here's the iterative version in Clojure:

(The N at the end of the numbers signifies that the numbers should have arbitrary precision: we might be dealing with some pretty big numbers here!)

Sure enough, we can generate the exact same set of 20 Fibonacci numbers that we had before.

     (map fibonacci (range 1 (inc 20)))
     ;> (1N 1N 2N 3N 5N 8N 13N 21N 34N 55N 89N 144N 233N 377N
        610N 987N 1597N 2584N 4181N 6765N)

This time though, we have a much more efficient implementation. Now we can actually calculate things like (fibonacci 500) with no trouble at all:

     (fibonacci 500)
     ;> 13942322456169788013972438287040728395007025658769730

        ..7264108962948325571622863290691557658876222521294125N

So let's analyze the time-complexity and space-complexity for this iterative process.

In this iterative method (shown by both the Python and Clojure versions above), note that we're looping up to n, and only doing some basic checks/variable assignment within that loop. So we have a process that grows only linearly with the input n, ie. Order of growth O(n) in time. 

And in terms of space, we only ever need space to hold the current and next values at any point! Now of course as the number of digits in the numbers increases we'll need more memory to store them, but in terms of algorithm analysis we consider this basically constant. So we have O(1) space-complexity.

In terms of algorithms, it doesn't get much better than linear time – if you've got a linear time algorithm you're doing pretty well, really. But can we do even better with the Fibonacci numbers? As it turns out, yes!

3. Using a Logarithmic-time Process

As noted in that good old SICP masterpiece, there is a way to compute the Fibonacci numbers with a logarithmic-time process. This algorithm is one of the coolest things I've ever seen, and it's also got a hefty bit of mathematics to it too!

So in the iterative process version that we had just before, we had the state transformation [current, next] → [next, current+next]. Let's denote this transformation by T. As we saw, applying T n times to the pair [0,1] gave us the nth Fibonacci number.

So mathematically, we could write that as Tⁿ(0,1) = [F_n , F_n+1].

Now here's the trick; suppose we consider T to be a special case of a more general transformation? Suppose that T_pq is a transformation that maps 
[current, next] → [p(current) + q(next), q(current) + (p+q)(next)].

  

Then you can see that taking p=0, q=1 gives us our Fibonacci transformation T.

Now consider what happens when we compose T_pq with itself:

So we have that T_pq² = T_p'q' , which gives us a notion of “squaring” this transformation that we can use to cut down on the amount of work the iterative process does! 

This allows us to speed up the iterative process we had earlier. We'll use the same loop-recur construct (modifying the values for the current, next parameters as we count down from n), but this time we'll check whether n is even at each point. If it is, then we can save some steps by applying the above "squaring" transformation instead. 

Here's some Clojure code that implements this new algorithm:

The code works exactly the same way as before, giving us the same output.

     (map fibonacci (range 1 (inc 20)))
     ;> (1N 1N 2N 3N 5N 8N 13N 21N 34N 55N 89N 144N 233N 377N 

        610N 987N 1597N 2584N 4181N 6765N)

However, this algorithm is efficient enough to calculate some huge Fibonacci numbers in no time at all, too. The following calculation took about half a second, and a good majority of that would've been printing it to the screen:

     (fibonacci 100000)
     ;> 25974069347...(absolutely huge number that spans
        over 10 A4 pages on a 12pt font)...3428746875N


We can modify the Python version of the iterative code to use this new logarithmic algorithm too:



The Python version works just as quick as the Clojure version:

     map(fibonacci, range(1,21))
     #> [1, 1, 2, 3, 5, 8, 13, 21, 34, 55, 89, 144, 233, 377, 
        610, 987, 1597, 2584, 4181, 6765]

     fibonacci(100000)
     #> 25974069347221...lots of digits...3428746875L

In terms of time- and space-complexity, we have a logarithmic order of growth in time (hence the name). Note too that we've still got that underlying iteration process happening – we're not delaying any operations or anything. So we're still working with constant space here. As such, we have O(log n) order of growth for time, and O(1) for space.

So seems like we've reached the epitome of efficiency then. We've gone from an exponential time process, to a linear time process, to a logarithmic time process. That's pretty awesome as it is, and the only thing more impressive would be if we could somehow compute a Fibonacci number in constant time, using some sort of mathematical formula. But there's no way that could be possible, right?

 

4. Using a Mathematical Formula

Remember that Golden Ratio we saw at the start? What if I said that, mathematically, this is true:

Of course, you shouldn't take my word for it! If we're going to suggest that this is true, and then go ahead and use it in code, we'll need to prove it first.

So here's a proof, (as per the one hinted at in the book):

 

So good for us! This means that we can use this formula to calculate the Fibonacci numbers directly, as we've done in the following Clojure function.

(Note that we've used Clojure's numeric-tower math library here to handle the rounding, exponentiation and square root stuff.)

Sure enough, this function works exactly like the others... right?

     (map fibonacci (range 1 (inc 20))) 
     ;> (1N 1N 2N 3N 5N 8N 13N 21N 34N 55N 89N 144N 233N 377N 

        610N 987N 1597N 2584N 4181N 6765N)

     (fibonacci 50) 
     ;> 12586269025N

     (fibonacci 70) 
     ;> 190392490709135N

     (fibonacci 90) 
     ;> 2880067194370824700N

     (fibonacci 100) 
     ;> 354224848179263100000N

     (fibonacci 200) 
     ;> 280571172992512000000000000000000000000000N

     (fibonacci 300) 
     ;> 222232244629422680000000000000000000000

        ..000000000000000000000000N

We're starting to see a suspicious amount of zeros at the end of these numbers...

     (fibonacci 100000)
     ;> NumberFormatException   java.math.BigDecimal.<init>
        (BigDecimal.java:470)

And for this one, we even get a cryptic error for our trouble!
A little step-by-step debugging reveals the culprit almost immediately:

     (math/expt golden-ratio 100000)
     ;> Infinity

So we pass in Infinity to the bigdec function. No wonder it complains!

The problem here is in the exponentiation, which uses Java Doubles. The second we pass the limit of that representation:

     (math/expt golden-ratio 1474)
     ;> 1.1163020658835235E308

     (math/expt golden-ratio 1475)
     ;> Infinity

     java.lang.Double/MAX_VALUE
     ;> 1.7976931348623157E308

.. we'll get Infinity from then on.

Also, using Doubles explains the zeros we were seeing earlier. As you can see from the scientific-notation form, they can only have so much precision.

Now we could attempt to solve both of those issues by using arbitrary precision decimals. However keep in mind that, compared with the other methods where we were only dealing with integers, this time we're dealing with floating-point approximations to real numbers. In particular, the operations that we're performing here – exponentiation, square-rooting, and even the simple division – will not work so well with arbitrary precision (after all, how could they? The whole point of many real numbers is that they don't have a terminating decimal expansion. For example, is the square root of 2 going to be 1.414214 ? Or 1.4142136 ? Or 1.41421356 ? Or ...)

But we actually have another problem too. This formula isn't the computational panacea that our earlier proof might have suggested. It's probably not possible to raise a number to the nth power in O(1) time, and judging by the algorithm used by the numeric-tower library we have at least
O(log n) time-complexity anyway.

That's still pretty awesome in terms of computational efficiency, but we already had a logarithmic-time algorithm before, didn't we? And in that one, we had arbitrary precision integers too!

This method is more interesting in a mathematical sense then, but not so practical for actual use.

So we've seen four different ways to generate the Fibonacci numbers here, and saw a few of the subtleties of algorithm analysis and space/memory concerns on our hunt towards the most efficient method we could find. We even saw some mathematics being put to good use here! It's all very interesting, isn't it?