Taylor series with Python and Sympy

Here I am again using my beloved Python and doing maths stuff.

Today I’d like to post a short piece of code I made after a review of Taylor series I did. This script lets you input (almost) any function, provided that it can be represented using Sympy and output the Taylor series of that function up to the nth term centred at x0.

Sympy is a great module for basic symbolic mathematics, it works fine and it is really simple to use even if you are new to Python.

Here is the output of the plot function for the function sin(x) approximating up to the 9th term:

	# >>>
	# Taylor expansion at n=1 x
	# Taylor expansion at n=3 -x**3/6 + x
	# Taylor expansion at n=5 x**5/120 - x**3/6 + x
	# Taylor expansion at n=7 -x**7/5040 + x**5/120 - x**3/6 + x
	# Taylor expansion at n=9 x**9/362880 - x**7/5040 + x**5/120 - x**3/6 + x

view raw taylor_series_out.py hosted with ❤ by GitHub

Note that in the console output the series is written backwards, however I think it could be fixed.

Here is the code I used

	import sympy as sy
	import numpy as np
	from sympy.functions import sin,cos
	import matplotlib.pyplot as plt
	plt.style.use("ggplot")
	# Define the variable and the function to approximate
	x = sy.Symbol('x')
	f = sin(x)
	# Factorial function
	def factorial(n):
	if n <= 0:
	return 1
	else:
	return n*factorial(n-1)
	# Taylor approximation at x0 of the function 'function'
	def taylor(function,x0,n):
	i = 0
	p = 0
	while i <= n:
	p = p + (function.diff(x,i).subs(x,x0))/(factorial(i))*(x-x0)**i
	i += 1
	return p

view raw taylor_series_p2.py hosted with ❤ by GitHub

This code can be customized to return Taylor expansions for any function you’d like to use (of course provided it can be represented using Sympy).

Finally, the code used to generate the plot

	# Plot results
	def plot():
	x_lims = [-5,5]
	x1 = np.linspace(x_lims[0],x_lims[1],800)
	y1 = []
	# Approximate up until 10 starting from 1 and using steps of 2
	for j in range(1,10,2):
	func = taylor(f,0,j)
	print('Taylor expansion at n='+str(j),func)
	for k in x1:
	y1.append(func.subs(x,k))
	plt.plot(x1,y1,label='order '+str(j))
	y1 = []
	# Plot the function to approximate (sine, in this case)
	plt.plot(x1,np.sin(x1),label='sin of x')
	plt.xlim(x_lims)
	plt.ylim([-5,5])
	plt.xlabel('x')
	plt.ylabel('y')
	plt.legend()
	plt.grid(True)
	plt.title('Taylor series approximation')
	plt.show()
	plot()

view raw taylor_series_p3.py hosted with ❤ by GitHub