Fibonacci Numbers - Golden Ratio (float)

In the previous post we discovered the relation between fibonacci numbers and golden ratio. This time, we attempt to calculate nth fibonacci number using the golden ratio.

/****************************************************************
 Program to generate nth Fibonacci number using Golden Ratio
 
 Theory: From Cormen
  fibonacci numbers are related to the golden ratio and 

  to its conjugate, which are the two roots of the

  equation
   x^2 = x + 1
  The roots are 
   r1=(1+sqrt(5))/2 = 1.61803... and
   r2=(1-sqrt(5))/2 = 0.61803...
   
  For nth Fibonacci number, we have
  F(n) = ((r1)^n - (r2)^n)/sqrt(5)
  
  Since r2<1, the value of r2 diminishes as the power
  applied on r2 increases. By induction, we can prove
  r2/sqrt(5) < 1/sqrt(5) < 1/2
  
  So finally we have 
   F(n) = floor(((r1)^n)/sqrt(5) + 0.5)
   
  Our approximation
   F(n) = floor((r1)^n)/sqrt(5)  

  works well in practice
   
 In our program we will use the same golden ratio to produce
 nth fibonacci number. 
 
 The limit when using the long data type for holding the 
 fibonacci result is 46th term.

 With float this becomes more interesting. 
 
***************************************************************/

#include <stdio.h>
#include <math.h>

const double sqrtOf5 = 2.2360679774997896964091736687313; //sqrt(5)
const double phi = 1.6180339887498948482045868343656; //(double)(sqrtOf5+1)/2; // we precalculate the root

int main(int argc, char* argv[])
{
 int n;      // out input value
 int i;
 
 double fibnum, fibfloor;
 double fibres;
 
 // get input from user
 printf("Program to print nth number of Fibonacci Series\n");
 printf("Enter n:");
 scanf("%d", &n);
 
 if (n <= 0)
 {
  printf("Invalid value of n. Please enter value > 0\n");
  goto exit;
 }

 printf("The nth Fibonacci Number is : ");
 
 fibnum = (pow(phi, n))/sqrtOf5;
 fibfloor = floor(fibnum);
 fibres = (fibnum - fibfloor)>0.5? ceil(fibnum):fibfloor; // approximate to nearest num
  
 printf("%.0f \n", fibres);
 
exit:
 printf("Exit from Program\n");
 
 return 0;
}

/***************************************************************
output:
E:\VSProjs\blog>2c_Fibonacci_nth_float.exe
Program to print nth number of Fibonacci Series
Enter n:76
The nth Fibonacci Number is : 3416454622906715
Exit from Program

E:\VSProjs\blog>2c_Fibonacci_nth_float.exe
Program to print nth number of Fibonacci Series
Enter n:77
The nth Fibonacci Number is : 5527939700884771
Exit from Program

E:\VSProjs\blog>2c_Fibonacci_nth_float.exe
Program to print nth number of Fibonacci Series
Enter n:78
The nth Fibonacci Number is : 8944394323791488
Exit from Program

E:\VSProjs\blog>2c_Fibonacci_nth_float.exe
Program to print nth number of Fibonacci Series
Enter n:79
The nth Fibonacci Number is : 14472334024676260
Exit from Program

PS E:\VSProjs\blog> .\2c_Fibonacci_nth_float.exe
Program to print nth number of Fibonacci Series
Enter n:82
The nth Fibonacci Number is : 61305790721611752
Exit from Program

PS E:\VSProjs\blog> .\2c_Fibonacci_nth_float.exe
Program to print nth number of Fibonacci Series
Enter n:83
The nth Fibonacci Number is : 99194853094755776
Exit from Program

PS E:\VSProjs\blog> .\2c_Fibonacci_nth_float.exe
Program to print nth number of Fibonacci Series
Enter n:84
The nth Fibonacci Number is : 160500643816367550
Exit from Program
***************************************************************/

Notes :

 We see that the fib(79) != fib(78) + fib (77). The expected
answer is : 14472334024676259. If we observe closely, the value
has been approximated (i.e the last two digits 59 is approxmiated
to 60. 

    But if we continue further, we get the right answer for 82
and 83, and again it flips at 84. 

 This is interesting. Everybody must already be aware that the 
way floating points are represented in computers are not the same as 
an integer. 
 
 So we have hit the limit for float. After this we lose precision
and zeros are attached. 

Representation of floating number
        The way floating numbers are represented in computers are by 
splitting the number into three parts. 
     - precision
     - exponent
     - sign
 (essentially as (sign)*(Mantissa or precision)*base^exponent
 eg. 13456*10^-2 represents 134.56)
 
 The IEEE 754 standard specifies a binary32(there is also a 
decimal32 format) as having:
  Sign bit: 1 bit
  Exponent width: 8 bits
  Significand precision: 24 (23 explicitly stored)
  
 s eeeeeeee fffffffffffffffffffffff
 
 With this format, we use 23 bits to store the number in binary
format. What we need to take a moment to realize is that the 
number is stored in binary format. That means the number of digits
of precision (when converted to decimal) vary. 
 
 eg. consider we have 4 bits to represent precision. we can 
represent numbers 0000 to 1111 in binary. Translated to decimal its
0 to 16. So the precision supported is 1 to 2. 
 
 Single precision number gives from 6 to 9 significant 
decimal digits precision. For double precision numbers it is
from 15 - 17 significant decimal digits. 
 
 The wikipedia article dedicated to this concept gives better
 insight into the format along with representing NaN, 
 Normalization, decimal32, binary32 etc.
 
Bottom Line is when using float or double, its all about precision 
and not the range of numbers you can represent. We can try and 
calculate fib(100). The value will not overflow like integer types, 
but we will not get the correct answer. 

For a problem where the outcome is dependent on the precision, we
need to be careful when choosing the floating number datatype.