Q10.1.4: Floating point approximations to two trigonometric functions
Question Q10.1.4
The functions $f(x) = (1-\cos^2x)/x^2$ and $g(x) = \sin^2x / x^2$ are mathematically identical, but plotted in the region $-0.001 \le x \le 0.001$ show a significant difference. Explain the origin of this difference.
Solution Q10.1.4
The expression 1 - np.cos(x)**2 suffers from catastrophic cancellation close to x=0 resulting in a dramatic loss of precision and wild oscillations in the plot of $f(x)$ (below). Consider, for example, x = 1.e-9: in this case, the difference1 - np.cos(x)**2 is indistinguishable from zero (at double precision) so f(x) returns 0. Conversely, np.sin(x)**2 is indistinguishable from x**2 and g(x) returns 1.0 correctly.
import numpy as np
import matplotlib.pyplot as plt
f = lambda x: (1 - np.cos(x) ** 2) / x**2
g = lambda x: (np.sin(x) / x) ** 2
x = np.linspace(-0.0001, 0.0001, 10000)
plt.plot(x, f(x), label=r"$(1-\cos^2x)/x^2$")
plt.plot(x, g(x), label=r"$\sin^2x/x^2$")
plt.ylim(0.99995, 1.00005)
plt.legend()
plt.show()