Examples for Week 5¶
Example 1:
It is known that for $N\gg{1}$,
$$
\dfrac{\pi^4}{90} \simeq \sum_{k=1}^N\dfrac{1}{k^4}\,.
$$
The larger $N$ is, the better the approximation of the value on the LHS by the sum on the RHS.
We want to find the relative errors of this approximation for $N=25, 50, 75, 100, 150, 200$
In [19]:
from math import *
def sum4(N):
sum = 1
for k in range(2,N):
sum += k**(-4)
return sum
def main():
exact = (pi**4)/90.0
Nvals = [25, 50, 75, 100, 150, 200]
for N in Nvals:
approx = sum4(N)
rel_error = fabs(exact - approx)/fabs(exact)
print("N={:<3d} {} R.E.= {:>6.2e}".format(N, "|", rel_error))
main()
Example 2: Use the central finite-difference formula to plot the derivative of the function $f(x) = \cos(x)$ over the interval $(0, \pi)$.
In [45]:
import math
import numpy as np
from matplotlib import pyplot as plt
# the function whose 1st derivative you want to approximate:
def my_f(x):
return np.cos(x) # for 'cos' you don't need this
# it is meant for more complicated functions
# second-order accuracy, 1st-order central finite diff. approx.
def ddf(fname,x,h=1.0E-5):
return (fname(x+h)-fname(x-h))/(2.0*h)
def main():
X = np.linspace(0, math.pi, 60)
N = len(X)
tmpY = []
for i in range(1, N-1):
tmp = ddf(my_f, X[i])
tmpY.append(tmp)
Y = np.array(tmpY)
X = X[1:N-1]
x = np.linspace(0,math.pi, 100)
y = -np.sin(x)
plt.plot(X,Y,'bo', label='CENTRAL FD')
plt.plot(x,y,'r', label='EXACT')
plt.xlabel('x')
plt.ylabel('y')
plt.legend(loc='upper center', fontsize='medium')
plt.title('Approximations for the 1st order derivative of COS(X)')
plt.show()
main()
Example 3: Use the Maclaurin approximation for the sine function to provide a means for calculating $\sin(x)$ to a prescribed accuracy, where $x$ is a user input.
The solution appears below:
In [20]:
"""
A simple implementation of the SINE function
"""
#
from math import fabs
import numpy as np
#--------------------------------------
def mySin(x):
# evaluates sin(x) from
# its (Maclaurin) power-series representation
EPS = 1e-11
i = 1
f, t = x, x
x2 = x*x
while (fabs(t) > EPS*fabs(f)):
i += 2 # this means i = i + 2
t *= -x2/float(i*(i-1.0))
f += t
return f
#--------------------------------------
# this is the demo function:
def main():
# test the function mySin on two cases:
# first case:
a1, a2 = np.sin(np.pi/4.0), mySin(np.pi/4.0)
rel_tol1 = fabs((a1-a2)/float(a2))
print('\n')
print('exact={:5.8f}, approximate={:5.8f}'.format(a1, a2))
print('\nRel. Error 1: {:5.5e}'.format(rel_tol1))
print('------------------------------------------------\n')
# second case:
b1, b2 = np.sin(1.5*np.pi), mySin(1.5*np.pi)
print('exact={:5.8f}, approximate={:5.8f}'.format(b1, b2))
rel_tol2 = fabs((b1-b2)/float(b2))
print ('\nRel. Error 2: {:5.5e}'.format(rel_tol2))
print('-------------------------------------------------\n')
#
# call the demo function:
main()