relation of the Fibonacci number series and the Golden ratio
Fibonacci number series:
Fibonacci number series displayed by a recurrence relation:
Galden Ratio:
https://en.wikipedia.org/wiki/Fibonacci_number
Fibonacci numbers are strongly related to the golden ratio: Binet's formula expresses the nth Fibonacci number in terms of n and the golden ratio, and implies that the ratio of two consecutive Fibonacci numbers tends to the golden ratio as n increases.
Fibonacci number series by Python
f = [0,1]
for _ in range(50):
    f+=[f[-1]+f[-2]]
print(f)
[0, 1, 1, 2, 3, 5, 8, 13, 21, 34, 55, 89, 144, 233, 377, 610, 987, 1597, 2584, 4181, 6765, 10946, 17711, 28657, 46368, 75025, 121393, 196418, 317811, 514229, 832040, 1346269, 2178309, 3524578, 5702887, 9227465, 14930352, 24157817, 39088169, 63245986, 102334155, 165580141, 267914296, 433494437, 701408733, 1134903170, 1836311903, 2971215073, 4807526976, 7778742049, 12586269025, 20365011074]
Find the Golen ratio from Fibonacci number series
Solve by SymPy
from sympy import Symbol,solve,factor
a = Symbol('A')
func = a*a + a - 1
a = factor(solve([func, a>0]))
A = float(a.rhs.n())
display(func,a,a.n())
A
0.6180339887498949
Numerical iteration
Get the golden ratio from the fibonacci
table = [(i,a/b) for i,(a,b) in enumerate(zip(f,f[1:]))]
table
[(0, 0.0),
 (1, 1.0),
 (2, 0.5),
 (3, 0.6666666666666666),
 (4, 0.6),
 (5, 0.625),
 (6, 0.6153846153846154),
 (7, 0.6190476190476191),
 (8, 0.6176470588235294),
 (9, 0.6181818181818182),
 (10, 0.6179775280898876),
 (11, 0.6180555555555556),
 (12, 0.6180257510729614),
 (13, 0.6180371352785146),
 (14, 0.6180327868852459),
 (15, 0.6180344478216818),
 (16, 0.6180338134001252),
 (17, 0.6180340557275542),
 (18, 0.6180339631667066),
 (19, 0.6180339985218034),
 (20, 0.618033985017358),
 (21, 0.6180339901755971),
 (22, 0.6180339882053251),
 (23, 0.618033988957902),
 (24, 0.6180339886704432),
 (25, 0.6180339887802427),
 (26, 0.618033988738303),
 (27, 0.6180339887543226),
 (28, 0.6180339887482036),
 (29, 0.6180339887505408),
 (30, 0.6180339887496481),
 (31, 0.618033988749989),
 (32, 0.6180339887498588),
 (33, 0.6180339887499086),
 (34, 0.6180339887498896),
 (35, 0.6180339887498969),
 (36, 0.6180339887498941),
 (37, 0.6180339887498951),
 (38, 0.6180339887498948),
 (39, 0.6180339887498949),
 (40, 0.6180339887498948),
 (41, 0.6180339887498949),
 (42, 0.6180339887498948),
 (43, 0.6180339887498949),
 (44, 0.6180339887498949),
 (45, 0.6180339887498949),
 (46, 0.6180339887498949),
 (47, 0.6180339887498949),
 (48, 0.6180339887498949),
 (49, 0.6180339887498949),
 (50, 0.6180339887498949)]
Good. Then, give a bias by -A to see the errors
import math
error = [(i,a/b - A) for i,(a,b) in enumerate(zip(f,f[1:]))]
error
[(0, -0.6180339887498949),
 (1, 0.3819660112501051),
 (2, -0.1180339887498949),
 (3, 0.04863267791677173),
 (4, -0.018033988749894925),
 (5, 0.0069660112501050975),
 (6, -0.0026493733652794837),
 (7, 0.0010136302977241662),
 (8, -0.00038692992636546464),
 (9, 0.00014782943192326314),
 (10, -5.6460660007306984e-05),
 (11, 2.15668056606777e-05),
 (12, -8.237676933475768e-06),
 (13, 3.1465286196574738e-06),
 (14, -1.201864649025275e-06),
 (15, 4.590717869179528e-07),
 (16, -1.7534976970434712e-07),
 (17, 6.697765930763211e-08),
 (18, -2.5583188345557062e-08),
 (19, 9.771908504596638e-09),
 (20, -3.732536946188247e-09),
 (21, 1.4257022229458016e-09),
 (22, -5.445698336714599e-10),
 (23, 2.080070560239733e-10),
 (24, -7.945166746736732e-11),
 (25, 3.034783535582619e-11),
 (26, -1.1591949622413722e-11),
 (27, 4.427680444507587e-12),
 (28, -1.6913137557139635e-12),
 (29, 6.459277557269161e-13),
 (30, -2.468025783741723e-13),
 (31, 9.414691248821327e-14),
 (32, -3.608224830031759e-14),
 (33, 1.3655743202889425e-14),
 (34, -5.329070518200751e-15),
 (35, 1.9984014443252818e-15),
 (36, -7.771561172376096e-16),
 (37, 2.220446049250313e-16),
 (38, -1.1102230246251565e-16),
 (39, 0.0),
 (40, -1.1102230246251565e-16),
 (41, 0.0),
 (42, -1.1102230246251565e-16),
 (43, 0.0),
 (44, 0.0),
 (45, 0.0),
 (46, 0.0),
 (47, 0.0),
 (48, 0.0),
 (49, 0.0),
 (50, 0.0)]
In the above, it indicates that the 39th iteration meets the golden ratio in 
To avoid the underflow, use decimal module.
import decimal
F = [0,1]
for _ in range(500):
    F+=[F[-1]+F[-2]]
decimal.getcontext().prec = 220
error2 = [(i,decimal.Decimal(a)/decimal.Decimal(b) - (-decimal.Decimal(1) + decimal.Decimal(5).sqrt())/decimal.Decimal(2))
          for i,(a,b) in enumerate(zip(F,F[1:]))]
error2[-50:]
[(451, Decimal('2.6586319293364935194863016834099E-189')),
 (452, Decimal('-1.0155070334308318486203693038262E-189')),
 (453, Decimal('3.878891709560020263748062280696E-190')),
 (454, Decimal('-1.481604794371742305040493803816E-190')),
 (455, Decimal('5.65922673555206651373419130762E-191')),
 (456, Decimal('-2.16163226293877649079763588460E-191')),
 (457, Decimal('8.2567005326426295865871634629E-192')),
 (458, Decimal('-3.1537789685401238517851315416E-192')),
 (459, Decimal('1.2046363729777419687682311629E-192')),
 (460, Decimal('-4.601301503931020545195619462E-193')),
 (461, Decimal('1.757540782015641947904546765E-193')),
 (462, Decimal('-6.71320842115905298518020825E-194')),
 (463, Decimal('2.56421744332073947649515719E-194')),
 (464, Decimal('-9.7944390880316544430526323E-195')),
 (465, Decimal('3.7411428308875685642063259E-195')),
 (466, Decimal('-1.4289894046310512495663445E-195')),
 (467, Decimal('5.458253830055851844927085E-196')),
 (468, Decimal('-2.084867443857043039117801E-196')),
 (469, Decimal('7.96348501515277272426327E-197')),
 (470, Decimal('-3.04178060688788778161170E-197')),
 (471, Decimal('1.16185680551089062057193E-197')),
 (472, Decimal('-4.4378980964478408010398E-198')),
 (473, Decimal('1.6951262342346161974012E-198')),
 (474, Decimal('-6.474806062560077911627E-199')),
 (475, Decimal('2.473155845334071760878E-199')),
 (476, Decimal('-9.44661473442137370999E-200')),
 (477, Decimal('3.60828574992340352128E-200')),
 (478, Decimal('-1.37824251534883685374E-200')),
 (479, Decimal('5.2644179612310704005E-201')),
 (480, Decimal('-2.0108287302048426632E-201')),
 (481, Decimal('7.680682293834575900E-202')),
 (482, Decimal('-2.933759579455301059E-202')),
 (483, Decimal('1.120596444531327286E-202')),
 (484, Decimal('-4.28029754138680789E-203')),
 (485, Decimal('1.63492817884715091E-203')),
 (486, Decimal('-6.2448699515464475E-204')),
 (487, Decimal('2.3853280661678342E-204')),
 (488, Decimal('-9.111142469570543E-205')),
 (489, Decimal('3.480146747033295E-205')),
 (490, Decimal('-1.329297771529334E-205')),
 (491, Decimal('5.07746567554716E-206')),
 (492, Decimal('-1.93941931134804E-206')),
 (493, Decimal('7.4079225849706E-207')),
 (494, Decimal('-2.8295746414305E-207')),
 (495, Decimal('1.0808013393219E-207')),
 (496, Decimal('-4.128293765343E-208')),
 (497, Decimal('1.576867902819E-208')),
 (498, Decimal('-6.02309943106E-209')),
 (499, Decimal('2.30061926507E-209')),
 (500, Decimal('-8.7875836406E-210'))]
Plot iteration
The red dashed line indicates the golden ratio (A = 0.6180339887498949).
The curve indicates the convergence to A in the oscillation condition
from matplotlib import pyplot as plt
import math
plt.figure(figsize=(15,4),dpi=80,facecolor='white')
plt.plot([a/b for a,b in zip(f,f[1:])])
plt.axhline(y=A,color='r',linestyle='dashed',alpha=0.4)
plt.title(r'$\frac{a_n}{a_{n+1}}$')
plt.yticks(ticks=[0,0.5,A,1],labels=['0.0','0.5','A','1.0'])
plt.grid()
plt.show()

See the errors
plt.figure(figsize=(16,15),dpi=60,facecolor='white')
for i,(start,end) in enumerate([(0,51),(0,7),(6,13),(20,27),(36,43),(42,49)],start=1):
    plt.subplot(3,2,i)
    plt.plot(*list(zip(*error[start:end])))
    plt.axhline(y=0,color='r',linestyle='dashed',alpha=0.4)
    plt.title(r'$\frac{a_n}{a_{n+1}}-A$, (Error['+str(start)+':'+str(end)+'])')
    plt.grid()
    plt.tight_layout()
plt.show()

Spread significant figures
list error issues underflow state, but list error2 works fine by the decimal.Decimal method.
plt.figure(figsize=(16,20),dpi=60,facecolor='white')
for i,(start,end) in enumerate([(0,501),(0,7),(6,13),(20,27),(36,43),(42,49),(194,201),(494,501)],start=1):
    plt.subplot(4,2,i)
    plt.plot(*list(zip(*error2[start:end])))
    plt.axhline(y=0,color='r',linestyle='dashed',alpha=0.4)
    plt.title(r'$\frac{a_n}{a_{n+1}}-A$, (Error2['+str(start)+':'+str(end)+'])')
    plt.grid()
    plt.tight_layout()
plt.show()

Recommended Posts