In [26]:
# this line is needed for the web notebook only
%matplotlib inline

CIVL4250 Numerical Methods in Engineering

Dr Dorival Pedroso

Code available on https://github.com/cpmech/CIVL4250py

Fibonacci numbers

See, e.g. Wikipedia

Starting with the three numbers: 0, 1, 1, the Fibonacci sequence can be constructed in such a way the new number is a sum of the previous two. Mathematically, if the new number is $F_k$, then:

$F_k = F_{k-2} + F_{k-1}$

These numbers can be easily generated in Python; however we will ignore the first zero number for the sake of convenience at the moment. First, let's create the initial list


In [27]:
x = [1, 1]

The number of Fibonnaci numbers, except the first three, will be


In [28]:
n = 21

The two previous numbers in x can be accessed as follows


In [29]:
last_but_one = x[-1] # also "penultimate", "second to last", "next to last"
last_but_two = x[-2] # also "antepenultimate", "two before the last"
print 'last but one =', last_but_one
print 'last but two =', last_but_two


last but one = 1
last but two = 1

Now, the sequence can be easily generated by appending the addition of the two last numbers, successively


In [30]:
print 'x (before) =', x
for i in range(n):
    x.append(x[-2] + x[-1])
print 'x (after) =', x


x (before) = [1, 1]
x (after) = [1, 1, 2, 3, 5, 8, 13, 21, 34, 55, 89, 144, 233, 377, 610, 987, 1597, 2584, 4181, 6765, 10946, 17711, 28657]

Golden ratio

By dividing each Fibonacci number by its predecessor, the results converge to a nice number, the golden ratio given by $\frac{1 + \sqrt{5}}{2}$.

In our code, the easiest way to verify this is to first convert the x list into a proper vector (array of real numbers). This can be accomplished by using array from numpy.

The array function gets x and a flag such as dtype=float as input. The flag tells array to convert integers to floats. The following code is then written


In [31]:
import numpy as np
x = np.array(x, dtype=float)
print 'x (vector) =', x


x (vector) = [  1.00000000e+00   1.00000000e+00   2.00000000e+00   3.00000000e+00
   5.00000000e+00   8.00000000e+00   1.30000000e+01   2.10000000e+01
   3.40000000e+01   5.50000000e+01   8.90000000e+01   1.44000000e+02
   2.33000000e+02   3.77000000e+02   6.10000000e+02   9.87000000e+02
   1.59700000e+03   2.58400000e+03   4.18100000e+03   6.76500000e+03
   1.09460000e+04   1.77110000e+04   2.86570000e+04]

All numbers after the first one in x and up to the end of the sequence can be recovered by using the slice notation. This is achieved as follows


In [32]:
x_upper = x[1:]
print 'x_upper =', x_upper


x_upper = [  1.00000000e+00   2.00000000e+00   3.00000000e+00   5.00000000e+00
   8.00000000e+00   1.30000000e+01   2.10000000e+01   3.40000000e+01
   5.50000000e+01   8.90000000e+01   1.44000000e+02   2.33000000e+02
   3.77000000e+02   6.10000000e+02   9.87000000e+02   1.59700000e+03
   2.58400000e+03   4.18100000e+03   6.76500000e+03   1.09460000e+04
   1.77110000e+04   2.86570000e+04]

Similarly, all numbers from the beginning of the sequence up to the last but one are obtained by means of


In [33]:
x_lower = x[:-1]
print 'x_lower =', x_lower


x_lower = [  1.00000000e+00   1.00000000e+00   2.00000000e+00   3.00000000e+00
   5.00000000e+00   8.00000000e+00   1.30000000e+01   2.10000000e+01
   3.40000000e+01   5.50000000e+01   8.90000000e+01   1.44000000e+02
   2.33000000e+02   3.77000000e+02   6.10000000e+02   9.87000000e+02
   1.59700000e+03   2.58400000e+03   4.18100000e+03   6.76500000e+03
   1.09460000e+04   1.77110000e+04]

Now, we can divide each component of one array by the corresponding one in the other array with the following code


In [34]:
g = x_upper / x_lower
print 'g =', g


g = [ 1.          2.          1.5         1.66666667  1.6         1.625
  1.61538462  1.61904762  1.61764706  1.61818182  1.61797753  1.61805556
  1.61802575  1.61803714  1.61803279  1.61803445  1.61803381  1.61803406
  1.61803396  1.618034    1.61803399  1.61803399]

These results can be plotted in order to visualise the convergence towards the golden ratio


In [35]:
import matplotlib.pyplot as plt
plt.plot(g);