#### Exercise A.1: Determine the limit of a sequence

#### a)

import numpy as np

def seq_a(N):
	a = np.zeros(N+1)
	for n in range(N+1):
		a[n] = (7+1/(n+1)) / (3-1/(n+1)**2) 
	return a

a = seq_a(100)
for n in range(len(a)):	
	print(f"a_{n} = {a[n]:7.5f}" )
	
"""
Terminal > python sequence_limits.py
a_0 = 4.00000
a_1 = 2.72727
a_2 = 2.53846
a_3 = 2.46809
.... some output skipped ....
a_97 = 2.33682
a_98 = 2.33678
a_99 = 2.33674
a_100 = 2.33671
"""

# Limit = 7/3 = 2.33333.....

#### b) 

def seq_D(N):
	D = np.zeros(N+1)
	for n in range(N+1):
		D[n] = np.sin(2**(-n)) / 2**(-n)
	return D

D = seq_D(10)
for n in range(len(D)):	
	print(f"a_{n} = {D[n]:7.5f}" )

"""
Terminal > python sequence_limits.py
.... output from part a) skipped ....
D_0 = 0.84147
D_1 = 0.95885
D_2 = 0.98962
D_3 = 0.99740
D_4 = 0.99935
D_5 = 0.99984
D_6 = 0.99996
D_7 = 0.99999
D_8 = 1.00000
D_9 = 1.00000
D_10 = 1.00000
"""

#### c) 

import matplotlib.pyplot as plt

def D(f, x, N):
	D = np.zeros(N+1)
	for n in range(N+1):
		h = 2**(-n)
		D[n] = (f(x+h)-f(x)) / h
	return D

# Let f(x)=sin(x), x=0, N=80
f = np.sin
x = 0
N = 80
Dn = D(f, x, N)
plt.figure()
plt.plot(np.arange(0,N+1), Dn, 'o')
plt.show()

#### d) 

# Let f(x)=sin(x), x=pi, N=80
f = np.sin
x = np.pi
N = 80
Dn = D(f, x, N)
plt.figure()
plt.plot(np.arange(0,N+1), Dn, 'o')
plt.show()

# Expected limit was cos(pi) = -1

#### e) 

# Plot numerator and denominator in Dn for last example
def D2(f, x, N):
	numer = np.zeros(N+1)
	denom = np.zeros(N+1)
	for n in range(N+1):
		h = 2**(-n)
		numer[n] = f(x+h)-f(x)
		denom[n] = h
	return numer, denom

f = np.sin	
x = np.pi
N = 80
numer, denom = D2(f, x, N)
for n,d in zip(numer, denom):
	print(f"{n:12.10f}  {d:12.10f}")
	if n==0:
		print("Numerator is zero")
	if d==0:
		print("Denominator is zero")
	
# Problem: numerator becomes zero