import numpy as np
%pylab inline
import matplotlib.pyplot as plt

#Number of iterations
Niter = 10000

#Double precision constants (Exact solution)
A_d = 1.0
B_d = -2/3.
#Signle precision constants
A_s = 1.0000000
B_s = -0.666667

#Solution to n-th term
pn = lambda A, B, n: A + B*n

#Arrays for storing the iterations
p_d = []
p_s = []
narray = range(Niter)
for n in narray:
    p_d.append( pn( A_d, B_d, n ) )
    p_s.append( pn( A_s, B_s, n ) )

#Converting to numpy arrays
p_d = np.array(p_d)
p_s = np.array(p_s)

#Relative error
error = p_d - p_s
plt.plot( narray, error, "-", color="blue" )
plt.xlabel("Number of iterations $n$", fontsize=14)
plt.ylabel("Error $p_n-\hat{p}_n$", fontsize=14)
plt.grid(True)

#New library!!! 
#Time: this library allows to access directly to the computer time. 
#      Useful when time calculations are required.
import time as tm

#Maximum number of iterations
narray = 10**np.arange(1,7,0.1)

#Time arrays
t_u = []  #User
t_p = []  #Python

#Iterations
for n in narray:
    #Generating random numbers
    N = np.random.random(n)
    #-------------------------------------------------------------------
    # MANUAL SUMMATION
    #-------------------------------------------------------------------
    #Starting time counter for user
    start = tm.clock()
    #Adding the numbers manually
    result = 0
    for i in xrange(int(n)):
        result += N[i]
    #Finishing time counter for user
    end = tm.clock()
    #Storing result
    t_u.append( end-start )
    
    #-------------------------------------------------------------------
    # PYTHON SUMMATION
    #-------------------------------------------------------------------
    #Starting time counter for user
    start = tm.clock()
    #Adding the numbers using python
    result = np.sum( N )
    #Finishing time counter for user
    end = tm.clock()
    #Storing result
    t_p.append( end-start )
    
#Ploting
plt.semilogx( narray, t_u, "-", color="red", linewidth=2, label="Manual" )
plt.semilogx( narray, t_p, "-", color="blue", linewidth=2, label="Python" )
plt.xlabel("Number of taken numbers $N$", fontsize=14)
plt.ylabel("Computing time [seconds]", fontsize=14)
plt.grid(True)

#Number of iterations
Niter = 100

#Double precision constants (Exact solution)
A_d = 1.0
B_d = 0.

#Solution to n-th term
pn = lambda A, B, n: A*(3.0)**-n + B*(3.0)**n

#Arrays for storing the iterations
p_s = [1.000000,0.333333]
p_d = [1.,1/3.]

narray = range(Niter)
for n in narray[2:]:
    p_s.append( 10/3.*p_s[n-1]-p_s[n-2] )
    p_d.append( pn( A_d, B_d, n ) )

#Converting to numpy arrays
p_d = np.array(p_d)
p_s = np.array(p_s)

#Relative error
error = p_d - p_s
plt.semilogy( narray, error, "-", color="blue" )
plt.xlabel("Number of iterations $n$", fontsize=14)
plt.ylabel("Error $p_n-\hat{p}_n$", fontsize=14)
plt.grid(True)

N = np.arange(1,1e3,0.1)
plt.loglog( N, N**2, linewidth=2, label="$N^2$" )
plt.loglog( N, N*np.log(N), linewidth=2, label="$N\ \log N$" )
plt.legend( loc="upper left" )
plt.grid(True)