Some libraries that we will use, numpy for various useful functions and arrays, matplotlib for plotting results



In [1]:

    
import numpy as np
import matplotlib
matplotlib.use('nbagg')
import matplotlib.pyplot as plt

Taylor Series

Taylor series will be used throuhout this notebook to derive numerical steppers, a brief explanation of the method follows:

A Taylor series is a representaion of a function as an infinite sum of (polynomial) terms calculated from the derivatives of the function at a fixed point.

The Taylor series of a function $f(x)$ around the point $a$ is given by:

$$f(a) + \dfrac{f'(a)}{1\!}(x-a) + \dfrac{f''(a)}{2!}(x-a)^2 + \dfrac{f'''(a)}{3!}(x-a)^3\dots $$

Or more concisely:

$$ f(x) \approx \sum_{n=0}^\infty \dfrac{f^{(n)}(a)}{n!}(x-a)^n $$

Where $f^{(n)}$ is the nth derivative of f.

This expression becomes useful when the series is truncated at some finite $n$ to obtain an approximation to a function of order $n$ with an error of order $n+1$.

$$ f(x) \approx \sum_{n=0}^m \left( \dfrac{f^{(n)}(a)}{n!}(x-a)^n \right) +\mathcal{O}(m+1) $$

When the higher order terms are neglected an approximation of the function is obtained.

$$ f(x) \approx \sum_{n=0}^m \left( \dfrac{f^{(n)}(a)}{n!}(x-a)^n \right) \\ $$

An example of this technique is shown below.

Consider the function $f(x)=e^x$

the derivatives of this function are trivially $e^x$

The taylor series of this function around the point a is therefore found using the above formula.

$$ f(x) \approx \sum_{n=0}^m \left( \dfrac{e^{(n)}(a)}{n!}(x-a)^n \right)\\ $$

Now let us consider the case where a=0

$$ f(x) \approx \sum_{n=0}^m \left( \dfrac{e^0}{n!}(x-0)^n \right)\\ f(x) \approx \sum_{n=0}^m \left( \dfrac{x^n}{n!} \right)\\ $$

An comparison of the taylor series of the exponential function with the exact function is shown below



In [2]:

    
def TaylorExponential(x, m=3):
    f=x*0
    if(m<1):
        print("Error: m must be a positive integer greater than 0")
        return
    for n in range(m):
       f+=np.power(x,n)/np.math.factorial(n)
    return(f)

for m in range(2,5):
    x = np.linspace(0,1,1000)
    f = TaylorExponential(x,m=m)
    plt.plot(x,f, label="m= {0}".format(m))
plt.plot(x,np.exp(x), label="Exact solution")
plt.xlabel("x")
plt.ylabel("f(x)")
plt.legend(loc="best")
plt.show()



In [ ]: