momom4urice
6/23/2017 - 7:51 AM

Finite difference derivative with Numdifftool (Fornberg method)

Finite difference derivative with Numdifftool (Fornberg method)

import numpy as np
import numdifftools.fornberg as ndf
import matplotlib.pyplot as plt
from random import random


def plot(outfile, x, fx=None, dfx=None, dfx_ref=None):
	plt.figure()
	if fx is not None:
		plt.plot(x, fx, label=r"$f(x)$", marker="*")
	if dfx is not None:
		plt.plot(x, dfx, label=r"$f'(x)$ (calculated)", marker="+")
	if dfx_ref is not None:
		plt.plot(x, dfx_ref, label=r"$f'(x)$ (reference)", marker="x")
	plt.legend()
	plt.xlabel("$x$")
	plt.savefig(outfile)
	plt.close()


def noisify(x, multiplier=1.):
	xmin, xmax = x.min(), x.max()
	deltax = xmax - xmin
	n = len(x)
	h = deltax / float(n-1)

	for i in range(len(x)):
		x[i] = x[i] * (1. + (2. * random() - 1.) / (float(n-1)) * multiplier)

	return x


# 1st order differentiation with even grid spacing
x = np.linspace(-1, 1, 21)
fx = np.tanh(x)
dfx = ndf.fd_derivative(fx, x, n=1)
dfx_ref = 1. - np.power(np.tanh(x), 2)
plot("even_tanh.png", x, fx, dfx, dfx_ref)

x = np.linspace(-np.pi, np.pi, 21)
fx = np.sin(x)
dfx = ndf.fd_derivative(fx, x, n=1)
dfx_ref = np.cos(x)
plot("even_sin.png", x, fx, dfx, dfx_ref)

x = np.linspace(0, 1, 21)
fx = np.exp(x)
dfx = ndf.fd_derivative(fx, x, n=1)
plot("even_exp.png", x, fx, dfx)

# 1st order differentiation with uneven (randomized) grid spacing
x = np.sort(noisify(np.linspace(-1, 1, 21), 5.))
fx = np.tanh(x)
dfx = ndf.fd_derivative(fx, x, n=1)
dfx_ref = 1. - np.power(np.tanh(x), 2)
plot("uneven_tanh.png", x, fx, dfx, dfx_ref)

x = np.sort(noisify(np.linspace(-np.pi, np.pi, 21), 5.))
fx = np.sin(x)
dfx = ndf.fd_derivative(fx, x, n=1)
dfx_ref = np.cos(x)
plot("uneven_sin.png", x, fx, dfx, dfx_ref)

x = np.sort(noisify(np.linspace(0, 1, 21), 5.))
fx = np.exp(x)
dfx = ndf.fd_derivative(fx, x, n=1)
plot("uneven_exp.png", x, fx, dfx)