""""
Ex. A.1 from "A primer on...
Compute various sequences and compare with 
the analytical results. 
"""


import numpy as np
import matplotlib.pyplot as plt

#a)
def seq_a(N):
    index_set = range(0,N+1,2)
    a = np.zeros(len(index_set))

    for i in range(len(a)):
        n = index_set[i]
        a[i] = (7 + 1.0/(n+1))/(3 - 1.0/(n+1)**2)

    return a

print('a):')
print(seq_a(100))
print(f'Analytical limit = {7/3}')


#b) 
def seq_b(N):
    index_set = range(N+1)
    d = np.zeros(len(index_set))

    for n in index_set:
        h = 2 ** (-n)
        d[n] = np.sin(h)/h
    
    return d

print('b):')
print(seq_b(50))

#c)
"""
Make a general function which uses a finite difference to 
estimate the derivative of a function f in a point x, 
and computes a sequence of how this approximate derivative 
depends on the discretization variable h (delta x).
"""

def D(f, x, N):
    index_set = range(N+1)
    d = np.zeros(len(index_set))

    for n in index_set:
        h = 2**(-n)
        d[n] = (f(x + h) - f(x))/ h
    
    return d


#apply the function to estimate the derivative 
#of sin(x) for x = 0:
d1 = D(np.sin, 0, 80)
#plt.plot(d1, 'ro')
#plt.show()


#d)
#Now, apply estimate the derivative 
#of sin(x) for x = pi:
d2 = D(np.sin, np.pi, 80)
plt.plot(d2,'go')
plt.show()

#e)
"""
The approximation in d) fails because of roundoff errors. 
The numerator becomes f(pi + h) - f(pi), and when h is
effictively zero this is f(pi)-f(pi) = 0. So even if the
denominator is also very small, the result is zero.
"""

"""
Terminal> python sequence_limits.py 
a):
[4.         2.53846154 2.43243243 2.39726027 2.38016529 2.37016575
 2.36363636 2.35905045 2.3556582  2.35304991 2.35098336 2.34930643
 2.34791889 2.34675206 2.34575734 2.34489938 2.34415187 2.34349483
 2.34291281 2.34239369 2.34192781 2.34150739 2.34112611 2.34077875
 2.34046098 2.34016919 2.33990031 2.33965175 2.3394213  2.33920705
 2.33900735 2.33882076 2.33864605 2.3384821  2.33832797 2.33818278
 2.33804579 2.33791632 2.33779377 2.3376776  2.33756732 2.3374625
 2.33736274 2.33726768 2.33717701 2.33709041 2.33700763 2.33692842
 2.33685255 2.33677981 2.33671002]
Analytical limit = 2.3333333333333335
b):
[0.84147098 0.95885108 0.98961584 0.99739787 0.99934909 0.99983725
 0.99995931 0.99998983 0.99999746 0.99999936 0.99999984 0.99999996
 0.99999999 1.         1.         1.         1.         1.
 1.         1.         1.         1.         1.         1.
 1.         1.         1.         1.         1.         1.
 1.         1.         1.         1.         1.         1.
 1.         1.         1.         1.         1.         1.
 1.         1.         1.         1.         1.         1.
 1.         1.         1.        ] 
"""