Finite-Differencing in NumPy

Recall from Difference approximation to first derivative, that a second-order accurate approximation to the first-derivative is:

\[\left . \frac{df}{dx} \right |_{x_0} = \frac{f(x_0 + \Delta x) - f(x_0 - \Delta x)}{2\Delta x} + \mathcal{O}(\Delta x^2)\]

We can compute the derivative of a function in NumPy that is sampled as an array by simply using views that are shifted by one element. In terms of array indices, \(i\), this would be

\[\frac{dx}{dx} \sim \frac{f(x_{i+1}) - f(x_{i-1}))}{x_{i+1} - x_{i-1}}\]

Note that for an array with \(N\) elements, we’ll have \(N-2\) derivatives, because we need an extra point on each end.

Here’s the code:

Listing 145 derivative.py

import numpy as np

# create our sampled data, f(x) = sin(x)
x = np.linspace(0.0, 2.0 * np.pi, 25)
f = np.sin(x)

# let's do it manually first
dfdx = np.zeros_like(f)
for i in range(1, len(f)-1):
    dfdx[i] = (f[i+1] - f[i-1]) / (x[i+1] - x[i-1])

# now construct (f[i + 1] - f[i-1]) / (x[i+1] - x[i-1])
# using views
dfdx2 = (f[2:] - f[:-2]) / (x[2:] - x[:-2])

# finally output, and compare to the exact result, cos(x)
for xc, fc, fc2 in zip(x[1:-1], dfdx[1:-1], dfdx2):
    print(f"{xc:7.4f}  {fc:7.4f}  {fc2:7.4f}  {np.cos(xc):7.4f}")

Notice:

• We have 2 fewer elements in our dfdx2 array than in our f array, because of how the views work (they leave our 2 elements either at the start or end)

• We compute the difference two ways: the first dfdx uses an explicit loop. This will be slower if our array is the same.