[TECH] A O(log(n)) time and O(1) space algorithm to find out Nth Fibonacci number.

We often come across Fibonacci numbers in several contexts, if you see Knuth's "Concrete Mathematics" in the chapter on "Generating functions" one interesting example illustrates the possible brick patterns to build a wall whose width is 'n' using bricks of dimension 2x1. Also see the UVA problem here

We know a simple dynamic programming based algorithm can lead us to a O(n) time algorithm with O(1). But we can do much better than that in fact we can do it in O(log(n)) time and O(1) space. We can rewrite the dynamic programming as matrix formulation as follows.

Fn Fn-1

=

1 1 1 0

Fn-1 Fn-2

This is really great we can get the closed form for F_n by repeated substitution. In fact the following is the closed form for F_n is as follows.

Fn Fn-1

=

1 1 1 0

n-1

1 0

So the problem of finding the n^th Fibonacci number is reduces to evaluating X^n-1 where X is a 2x2 matrix as shown above. Now the question is how we can evaluate X^n-1 as efficiently as possible. In fact to compute any a^b we need exactly |_{_}log(n)_{_}| multiplications by making use of the bit representation of b = (2^i₁ + 2^i₂ ...) where i_j are the positions of bits which have 1 in the binary representation of b. I used this trick to get a O(log(n)) time and O(1) time algorithm to compute the n^th Fibonacci number. Please see the code below get this here

  1 /*Find the n^th fibonacci number in O(log(n)) 
  2  *time, using the following idea of matrix formulation 
  3  * 
  4  * [   F_n  ]    [ 1      1  ] [F_{n-1}]
  5  * [        ] = [         ] [       ]
  6  * [ F_{n-1}]    [ 1      0  ] [F_{n-2}]
  7  *
  8  * vamsik(at)engr(dot)uconn July 11, 2008
  9  */
 10 #include<stdio.h>
 11 #include<stdlib.h>
 12 /*Fib returns the n^th Fibonacci number*/
 13 unsigned int Fib(unsigned int n){
 14     unsigned int deg=1;
 15     unsigned int fib[2][2] = {{1,1},{1,0}};
 16     unsigned int fibAns[2][2];
 17     unsigned char initAns=0;
 18     unsigned int mul[4];
 19     if(n<=2){ return (n<2)?n:2;}
 20     while(deg < n-1){
 21          mul[0] = fib[0][0]*fib[0][0] + fib[0][1]*fib[1][0];
 22          mul[1] = fib[0][0]*fib[0][1] + fib[0][1]*fib[1][1];
 23          mul[2] = fib[1][0]*fib[0][0] + fib[1][1]*fib[1][0];
 24          mul[3] = fib[1][0]*fib[0][1] + fib[1][1]*fib[1][1];
 25          fib[0][0] = mul[0]; fib[0][1] = mul[1]; 
 26          fib[1][0] = mul[2]; fib[1][1] = mul[3];
 27          deg <<=1;
 28          if(deg & n-1){ /*i^th bit set*/
 29              if(!initAns){
 30                  fibAns[0][0] = fib[0][0]; fibAns[0][1] = fib[0][1]; 
 31                  fibAns[1][0] = fib[1][0]; fibAns[1][1] = fib[1][1];
 32                 initAns=1;
 33              }else{
 34                  mul[0] = fibAns[0][0]*fib[0][0] + fibAns[0][1]*fib[1][0];
 35                  mul[1] = fibAns[0][0]*fib[0][1] + fibAns[0][1]*fib[1][1];
 36                  mul[2] = fibAns[1][0]*fib[0][0] + fibAns[1][1]*fib[1][0];
 37                  mul[3] = fibAns[1][0]*fib[0][1] + fibAns[1][1]*fib[1][1];
 38                  fibAns[0][0] = mul[0]; fibAns[0][1] = mul[1]; 
 39                  fibAns[1][0] = mul[2]; fibAns[1][1] = mul[3];
 40              }
 41          }
 42     }
 43     if(n-1 & 1){
 44         return fibAns[0][0] + fibAns[0][1];
 45     }
 46     return fibAns[0][0];
 47 }
 48 int main(int argc,char **argv){
 49     unsigned int n;
 50     while(1){
 51         scanf("%u",&n);
 52         if(!n){
 53             exit(0);
 54         }
 55         printf("%u\n",n,Fib(n));
 56     }
 57 }
 58